Bipartite Fluctuations as a Probe of Many-Body Entanglement 
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We investigate in detail the behavior of the bipartite fluctuations of particle number N and spin 
5"^ in many-body quantum systems, focusing on systems where such U(l) charges are both conserved 
and fluctuate within subsystems due to exchange of charges between subsystems. We propose that 
the bipartite fluctuations are an effective tool for studying many-body physics, particularly its en- 
tanglement properties, in the same way that noise and Full Counting Statistics have been used in 
mesoscopic transport and cold atomic gases. For systems that can be mapped to a problem of 
non-interacting fermions we show that the fluctuations and higher-order cumulants fully encode the 
information needed to determine the entanglement entropy as well as the full entanglement spec- 
trum through the Renyi entropies. In this connection we derive a simple formula that explicitly 
relates the eigenvalues of the reduced density matrix to the Renyi entropies of integer order for any 
finite density matrix. In other systems, particularly in one dimension, the fluctuations are in many 
ways similar but not equivalent to the entanglement entropy. Fluctuations are tractable analytically, 
computable numerically in both density matrix renormalization group and quantum Monte Carlo 
calculations, and in principle accessible in condensed matter and cold atom experiments. In the con- 
text of quantum point contacts, measurement of the second charge cumulant showing a logarithmic 
dependence on time would constitute a strong indication of many-body entanglement. 

PACS numbers: 03.67.Mn, 05.30.-d, 05.70.Jlc, 71.10.Pm 



I. INTRODUCTION 

The study of quantum many-body systems has 
traditionally involved the analysis of ground-state 
and excited-state energies, correlation functions, and 
symmetry-breaking order parameters to characterize dif- 
ferent states of matter [1, 2]. Recently, however, there has 
been great interest in using a different feature to under- 
stand the properties of a quantum many-body system, 
namely entanglement. Entanglement is believed to be 
particularly relevant at zero temperature, especially for 
topologically ordered states for which conventional order 
parameters are not sufficient to fully define the state of 
the system [3-5]. 

The degree of entanglement in a system is often quan- 
tified using the bipartite entanglement entropy between 
a subsystem and the remainder of the system, most com- 
monly using the scaling of the entanglement entropy 
with subsystem size. Specifically, given a wave function 
1^') — usually the ground state but often a time-evolved 
wave function — the entanglement entropy is defined as 
the von Neuman entropy of the reduced density matrix 
PA = Ttb I^*) (^fj of subsystem A obtained by tracing out 
the degrees of freedom in the remainder of the system B. 
Thus 

S{pA) = -TT{pA\npA). (1.1) 
More generally, the a-Renyi entropies are defined as 

'So(/3A) = ^ln[Tr(pl)] (1.2) 




FIG. 1: (color online). Quantum many-body systems. We 
can divide the system into two parts A and the remainder 
B, where A can be (a) a single connected region or, more 
generally, (b) multiple disconnected regions. The result is 
the formation of a boundary, across which two subsystems 
of a local Hamiltonian interact, for example by exchanging 
particles (arrows). 



and reduce to the von Neumann entropy in the limit 
a — > 1. For a pure state Sa{pA) = Sa{pB) and it makes 
no difference which subsystem is analyzed; indeed, this 
symmetry is an important feature of entanglement en- 
tropy. The scaling of entanglement entropy with sub- 
system size has yielded a wealth of interesting results for 
both gapped and critical systems [3, 6-19]. Alternatively, 
detailed information about the structure of entanglement 
can be obtained by studying the full spectrum of eigen- 
values of the reduced density matrix, specifically through 
the eigenvalues of the "entanglement Hamiltonian" Ha 



defined by pA = e'^-^ /Tiie'"^) [20-32]. 

A distinguishing feature of using entanglement entropy 
to study many-body systems is the introduction of a 
boundary between the two subsystems A and B as in 
Fig. 1, which represents a change in perspective from 
the study of point-point correlators and local order pa- 
rameters. The special role played by the area of this 
boundary has long been recognized, and indeed was an 
early motivation for the study of entanglement entropy 
[33, 34], because most Hamiltonians of interest contain 
local interactions which in turn lead to correlations local- 
ized along the entire boundary. In many situations this 
is formalized in area laws for the entanglement entropy 
[35], while important violations have been studied since 
then [8, 12]. The important role played by the bound- 
ary, moreover, is not unique to entanglement entropy — all 
proposed measures of many-body entanglement are based 
on it. In systems where particles are exchanged between 
parts of the system, for example, all particle exchange 
occurs through the boundary and therefore one expects, 
and can confirm, that particle number fluctuations will 
similarly be sensitive to the area. Indeed, as we have 
shown previously [36] and further elaborate on in this 
work, for non-interacting systems of fermions in which 
entanglement originates precisely from the uncertainty 
of the position of particles, particle number fluctuations 
are equivalent to the entanglement entropy. 

Fundamentally, entanglement is what distinguishes 
quantum systems from classical ones, and the belief that 
studying the entanglement properties of a quantum sys- 
tem, whether by analyzing the scaling of the entangle- 
ment entropy with subsystem size or by investigating the 
more detailed entanglement spectrum, has led to many 
striking and beautiful insights. In particular, the recent 
interest in many-body entanglement was in large part 
motivated by the work in Refs. [33, 34] showing that the 
entanglement entropy of a subsystem in a quantum held 
theory is generally expected to scale as the area of the 
subsystem, and the work of Refs. [7, 37] showing that 
in conformally invariant one-dimensional gapless systems 
the entanglement entropy of an interacting many-body 
system exhibits universal logarithmic scaling with sub- 
system size, with the prefactor controlled by the central 
charge of the underlying conformal fleld theory (CFT). 
These results have been conflrmed by many numerical, 
and sometimes analytical, calculations for microscopic 
Hamiltonians, and indeed the question of which problems 
are tractable by classical simulation methods has itself 
been an important part of the investigation of many-body 
entanglement [38]. 

Except for systems that can be mapped to non- 
interacting fermions [36] as suggested in Refs. [39, 40] 
and discussed extensively in the present work, however, 
it has not been clear how entanglement entropy can be 
measured experimentally (similar ideas are discussed in 
Ref. [41]; a method for determining the lower bound is 
proposed in Ref. [35]; see Ref. [42] for a completely dif- 
ferent, promising approach). Even theoretically, entan- 



glement entropy is difficult to calculate analytically and 
in dimensions greater than one there are no good numer- 
ical methods for obtaining this quantity due to the fact 
that the entanglement entropy is a non-trivial function of 
the reduced density matrix and hence of the full ground- 
state wave function (however, some signiflcant progress 
has been made recently for the Renyi entropies in quan- 
tum Monte Carlo (QMC) [43-46]). This has prompted 
a search for other measures of entanglement that behave 
similarly to the entanglement entropy yet are easier to 
determine, primarily in QMC calculations as they can be 
applied to a wide variety of systems in any dimension. 

In this work we show that using the same setup shown 
in Fig. 1 but measuring instead the fluctuations of a con- 
served U(l) charge such as particle number A'^ and spin 
reveals important features of quantum many-body 
systems including their entanglement properties, quite 
similarly to the way in which Full Counting Statistics 
(FCS), the study of charge transfer across mesoscopic 
conductors, has been intensely analyzed in mesoscopic 
transport [47-50] and in cold atom systems [51, 52]. In- 
deed, our result for non-interacting fermions unequivo- 
cally demonstrates the importance of studying the full set 
of cumulants, beyond the fluctuations (noise) encoded in 
the second cumulant. The present work, moreover, may 
be understood as applying some of the same principles 
underlying the fleld of FCS to the study of the ground 
states of many-body Hamiltonians, with the boundary 
introduced by the bipartition playing the role of the scat- 
tering or interaction region (e.g., quantum dots) in meso- 
scopic systems. 

To be precise, suppose the total particle number is con- 
served and that the pure state |^') is an eigenstate of 
the particle number operator iVtot, -Atot 1^*) = Atot 1^") 
where Atot (without the caret) is a positive, real inte- 
ger. Since we will almost always be concerned with pure 
states at zero temperature in this work, we only consider 
pure states here. We then deflne the fluctuations as [53] 

T={{NA-{NA)r) (1.3) 

where the expectation value is taken with respect to 1 5') . 
We can generalize the fluctuations to the cumulants de- 
fined by 

Cn = {-idxrinxW\x=o, (1.4) 

where for the ground state of a many-body Hamiltonian 
we might define 

X(A) = (exp(zA7VA)) (1.5) 

as the generating function for particle number fluctua- 
tions within subsystem A, with J- = C2- In the context 
of two-terminal mesoscopic transport the time-dependent 
generating function is deflned as [47] 

x(A,t) = ^P„,(t)e^^'", (1.6) 
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where Pm{i) is the probabihty that m charges are trans- 
ferred from the source lead to the drain lead. As we will 
see below, in the framework of Fig. 1 the even cumulants 
are the appropriate generalization of the fluctuations be- 
cause of their invariance under exchange of subsystems at 
zero temperature. We note that the fluctuations (and cu- 
mulants) can always be expressed as a sum over the sub- 
system of the appropriate correlation functions of density 
operators [cf. Eq. (3.23)], which shows how the fluctu- 
ations may be considered as the quantity that bridges 
local correlation functions with the highly non-local and 
non-linear entanglement entropy. As such the fluctua- 
tions are expected to share some, but not all, essential 
features with the entanglement entropy. 

The fluctuations at zero temperature (and more gen- 
erally all even cumulants) share with the entanglement 
entropy two important properties: 

1. The fluctuations (and all cumulants) are zero for a 
product state \'^) = \^ a) ® For the product 
state, \^ a) and are separately eigenstates of 
Na and Ng, respectively, which can be seen by 
noting that if either \^ a) or is a superposition 
of states with different numbers of particles then so 
must the total state l^*) since N is additive. Thus 
there are no particle number fluctuations in A or 
B. 

Unlike the entanglement entropy, the converse of 
this statement is not always true: there exist states 
with zero fluctuations which are not product states. 
This can occur because one can construct a Hamil- 
tonian where A'^ and Nb are separately conserved 
from the start. Then all eigenstates of the Hamil- 
tonian will have zero fluctuations even though the 
entanglement entropy is nonzero. We note that 
if the fluctuations, i.e., the second cumulant, are 
nonzero, then all cumulants are zero and implies 
that Na = {Na) and Nb = {Nb)- However, unhke 
total particle number conservation, particle num- 
ber conservation within subsystems is rare (essen- 
tially, only density-density interactions are allowed 
between A and B as in Ref. [54]) and we can gen- 
erally expect the converse of this statement to hold 
in realistic systems. 

2. At zero temperature the fluctuations are symmet- 
ric between A and the remainder of the system 
B. Since Na + Nb = Ntot we have Na - {Na) = 
-{Nb - {Nb)) and therefore 

J'a = J'b- (i.r) 

More generally, switching A and B is equivalent to 
changing A — > — A in the generating function (1.5). 
Thus the even cumulants are invariant under this 
switch and C„ for even n are symmetric in A and 
B as well. 

This property is crucial to understanding why the 
fluctuations are generally expected to behave sim- 



ilarly to the entanglement entropy, in the follow- 
ing sense. In some of the earliest works on the 
scaling of the entanglement entropy with the area 
of the subsystem, it was suggested [34] that this 
could be understood from the fact that the entan- 
glement entropy of a pure state is symmetric be- 
tween the subsystem and the remainder of the sys- 
tem, which implies that the relevant correlations lie 
at the boundary shared by the two systems and not 
in the "volume" of the respective subsystems. We 
now know from many violations of the area law that 
symmetry is not sufficient to guarantee an area law, 
but symmetry is a property shared by all proposed 
entanglement measures. 

The fluctuations = C2, in addition, share with the 
entanglement entropy the property of being subadditive, 
in the sense that (at zero temperature) 

J'a + J'b> 0, (1.8) 

which must be true because J-a = J'b > 0. The higher- 
order cumulants and Renyi entropies for a > 1 are not 
subadditive. 

Thus with the exception noted under Property 1 the 
fluctuations can, in a broad sense, be considered an en- 
tanglement measure similar to previously proposed mea- 
sures (disussed below), while the question of how simi- 
larly the fluctuations behave to the entanglement entropy 
itself is a separate and interesting question addressed in 
the following sections. It is, however, important to note 
that the fluctuations are not an entanglement measure 
in the strict quantum information sense, which is a non- 
negative function of a state which cannot increase under 
local operations and classical communications and which 
vanishes for separable states. On the other hand, given 
the prior knowledge that the ground state is an eigenstate 
of the particle number operator, for example, the fluctu- 
ations contain (perhaps incomplete) information about 
entanglement in the system. 

The advantage of the fluctuations, meanwhile, is that 
they are much more accessible — analytically, numeri- 
cally, and experimentally. Computing the fluctuations 
involves trivial adjustments to existing density matrix 
renormalization group (DMRG) [38, 55] and QMC codes 
[in DMRG the reduced density matrix is almost always 
stored in block-diagonal form according to the U(l) quan- 
tum number, so that the sum of the eigenvalues within 
each block can be used to easily determine cumulants of 
any order]. There are, moreover, existing experimental 
techniques that can in principle be used to measure the 
fluctuations in real systems [56-60] . This is the main mo- 
tivation of the present work: to explore the behavior of 
the fluctuations in a wide variety of systems in anticipa- 
tion of such experimental measurements. 

The outline of this paper is as follows. In Sec. II we 
show that for systems that can be mapped to a problem 
of non-interacting fermions the full set of even cumu- 
lants is completely equivalent to the entanglement en- 
tropy by showing how the von Neumann entanglement 
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entropy [36], and more generally the Renyi entanglement 
entropies, can be written exactly in terms of the even 
cumulants. In particular, this usually implies that the 
dominant scaling of any given cumulant with subsystem 
size is identical to that of the entanglement entropy. We 
present applications of the formulas, as well as a (com- 
pletely general) method for obtaining the full entangle- 
ment spectrum from the Renyi entropies, or in the case of 
non-interacting fermions, from the cumulants. In Sec. Ill 
we explore the generalization for one-dimensional inter- 
acting systems as suggested in Ref. [53], including the 
case of disjoint intervals (Fig. lb). Although the ex- 
act relation found in systems of non-interacting fermions 
no longer holds we show how many similarities exist be- 
tween fluctuations and the entanglement entropy. We 
then establish that for gapped quantum many-body sys- 
tems an area law holds for the fluctuations. Finally, in 
Sec. IV we compare bipartite fluctuations to previously 
proposed measures of many-body entanglement and show 
that fluctuations possess almost all of the same proper- 
ties while being much easier to access both theoretically 
and experimentally. We conclude by suggesting possible 
experimental realizations and future work. 

Appendices provide essential but somewhat technical 
calculations, while the Supplementary Material contains 
longer, more detailed computations for several of the re- 
sults presented in the main text. 

Unless stated otherwise, we work in units where the 
Planck constant h and Boltzmann constant ks are h = 
fcs = 1. For lattice systems we assume that the lattice 
spacing a = 1, while for continuum systems we will ex- 
plicitly include the short-distance cutoff. 



II. EXACT RELATIONS BETWEEN CHARGE 
STATISTICS AND ENTANGLEMENT 
ENTROPIES 



The restriction to non-interacting fermions is not as 
limiting as it flrst appears; aside from the intrinsic 
importance of free fermions in various dimensions, in 
one dimension the Jordan- Wigner transformation maps 
the spin- 1/2 XY model to a model of non-interacting 
fermions, and a similar mapping can be carried out for 
bosons with infinite on-site repulsion (hard-core bosons) 
to a model of free fermions. In two dimensions the in- 
teger quantum Hall effect is described by a model of 
non-interacting fermions, as are the related topological 
insulators. The case of hard-core bosons in one dimen- 
sion seems particularly promising for the experimental 
detection of many-body entanglement, since the num- 
ber of particles in a given region is already accessible in 
present-day experiments with cold atoms in optical lat- 
tices [61-63]. 

Before going through the detailed derivations, we flrst 
present the central results of this section. As we demon- 
strate below, the von Neumann entanglement entropy S 
can be expressed in terms of the charge statistics as [36] 



A'+l 



S = lim 



J2 (^n{K)C„ 



(2.1) 



where C„ are the cumulants deflned in Eqs. (1.4), (1.5). 
The cutoff-dependent coefficients an{K) are given by 



a„{K) 



^ Z^fc=n-l 





Sl(fc,71-1) 

k'.k 



for 71 even, 
for 71 odd, 



(2.2) 



where Si{n,m) are the unsigned Stirling numbers of the 
first kind. 

A similar, but lengthier, relation also exists for the 
Renyi entanglement entropies of integer order. 



lim 



J2l3kin,R)Ck. 



(2.3) 



In this section we show how, for non-interacting 
fermions (and systems that can be mapped to a problem 
of non-interacting fermions), two important measures of 
many-body entanglement can be expressed exactly in 
terms of the statistics of charge fluctuations. Speciflcally, 
we derive series for the von Neumann and the Renyi en- 
tanglement entropies in terms of the full set of (even) 
cumulants of the charge fluctuations. Similar ideas have 
previously been suggested in Ref. [39] in the context of 
ground-state entanglement entropy and in Ref. [40] in 
relation to electron transport through a quantum point 
contact (QPC). Here, we extend the results of Ref. [36] to 
the computation of the generalized Renyi entropies and 
demonstrate that such relations in fact hold in general 
for non-interacting fermionic systems in any dimension. 
Importantly, the intimate connections between the en- 
tanglement entropies and the charge fluctuations facili- 
tate a novel approach to the experimental determination 
of many-body entanglement entropy in non-interacting 
systems where charge fluctuations may be measured. 



Here, the coefficients depend on the cutoff R as 



f3kin,R) = < 



B. 



r, 

r \ r / \in 



-1 



Enr 
s=k 
/ nm \ Si(s,k) 
\7ir—s/ (s — 1)! 



r+s+Tir-f-nm 



for k even, 







for A: odd. 

(2.4) 

Interestingly, only even-order cumulants contribute in 
the series (2.1) and (2.3). Physically, this is a con- 
sequence of the entanglement entropies being symmet- 
ric between the subsystem and its complement for pure 
states. One may for example think of a QPC at zero tem- 
perature, Eq. (1.6), where Pm describes the probability 
of transmitting tti charges from one subsystem, the left 
lead, to its complement, the right lead: charge conser- 
vation implies that the number of electrons tti collected 
in the right lead is equal to —m charges collected in the 
left. Only even-order cumulants are then symmetric in 
the two leads, and the entanglement entropies for the 
QPC can therefore only depend on the even cumulants. 
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FIG. 2: (color online). Entanglement entropy of a pair of 
spins. The two spins are in the state \D) in Eq. (2.5), which is 
a product state for D = 0,1 and maximally entangled for D = 
1/2. The red line shows the exact result given by Eq. (2.6), 
while the blue lines show the series (2.1) with increasing cutoff 
number K from bottom to top. For K = 10, the series has 
basically converged to the exact result. 



Since the series (2.1) and (2.3) are universally valid for 
any non-interacting system, the entanglement entropies 
can only depend on the even-order cumulants, including 
for systems where the even cumulants are not necessarily 
symmetric under interchange of the two subsystems. We 
consider an example of such a system in Sec. II C 2. 

The two series above yield a sequence of approxima- 
tions to the entanglement entropies that is an increas- 
ingly better lower bound to the exact results: if the se- 
ries (2.1) for the von Neumann entanglement entropy is 
truncated at a finite K, the approximation only becomes 
better as K is increased. Thus, with only a finite num- 
ber of cumulants, the series provides a converging ap- 
proximation of the exact von Neumann entropy. This 
is also true for the Renyi entropies. For infinite cutoff 
K, our series (2.1) is moreover formally equivalent to 
the series suggested in Ref . [40] , which can be written as 
S = '}Z'^=i o-nCn, where q„ = 2^(n) for n even with C(n) 
being the Riemann zeta function, and zero for n odd. 
This is shown explicitly in the Supplementary Material. 
As we show below, however, the convergence properties 
of the series (2.1) and the one presented in Ref. [40] are 
very different. 

As a simple illustration of how the series can be ap- 
plied, consider two electron spins in the state 



\D)ab = yr^l t)A| i)B - 1)a\ t)j 



(2.5) 



with < D < 1. Clearly, for D = 1/2 the spins are in 
the maximally entangled singlet state, while for D = 0, 1 
the spin state is factorized. It is easy to calculate the 
entanglement entropy of spin A with respect to spin B 
(or vice versa) directly from the definition (1.1), 

S = -[DlnD + {l-D)ln{l-D)], (2.6) 

which is zero for D = 0.1 and the maximal In 2 for 



D = 1/2. The entanglement entropy can also be ob- 
tained from the fiuctuations of the spins using the series 
(2.1). For D = 0, 1 the two spins are clearly in a well- 
defined state and no fluctuations occur. In contrast, in 
the spin singlet state D = 1/2 the two spin states of spin 
A are equally probable and the 2-component of each spin 
fluctuates. We can think of the two spins as sharing a 
spin excitation: for the state ] t)yi| 4-)b spin A is excited, 
while for the state | t) b the excitation has been 

transferred to spin B. 

Now, the distribution P(m) for the number m of trans- 
ferred excitations between spin A and spin B is 



P(m) 



l-D for m = 0, 
D for m = 1, 

otherwise. 



The corresponding generating function is 

X(A) = V P(m)e'^™ = 1 + (e^^ - l)D, 



(2.7) 



(2.8) 



from which the cumulants of m easily follow as C„ = 
(— ic^a)" lnx(A)|A=o- In Fig. 2 we show the entanglement 
entropy calculated from the cumulants using the series 
(2.1) with increasing cutoff number K together with the 
exact result, Eq. (2.6). Already with K = 2, where only 
the second cumulant is included, the series yields a good 
approximation of the exact entanglement entropy, while 
basically full convergence is obtained for K = 10. We 
note that the convergence is not uniform over all values 
of D. This simple example illustrates how the entan- 
glement entropy can be obtained from the (correlated) 
fluctuations of the spins. 

It is instructive to attempt the same calculation us- 
ing the series presented in Ref. [40]. The terms of this 
series are 2(^{n)C„, where 2({n) ~ 2 for large n. How- 
ever, the cumulants corresponding to Eq. (2.8) go as 
\Cn\ ~ {"IT- ~ 1) ^01' large n, so that the series of 
Ref. [40] diverges under normal summation methods. 
Since the factorial divergence of the |C„| is typical for 
non-gaussian fluctuations [60], moreover, the problem is 
quite generic. As explained in the Supplementary Ma- 
terial, this is due to the existence of singularities in the 
complex plane possessed by almost all cumulant gener- 
ating functions and hence the problematic step of inte- 
grating its series expansion term by term. In particular, 
for any experimental situation where only a finite num- 
ber of cumulants is available, the series from Ref. [40] 
is ill-suited for the determination of the von Neumann 
entanglement entropy. Instead the series (2.1) should be 
used as it converges for any distribution of the fluctua- 
tions. 



A. Derivation 

We now present the detailed derivation of the se- 
ries (2.1) for the von Neumann entanglement entropy. 
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The series (2.3) for the Renyi entropies follows from a 
similar but lengthier derivation, an outline of which is 
given in Appendix A. In Sec. II C we illustrate the use 
of the series with several different applications, and the 
reader may at this point wish to proceed directly to the 
concrete examples which can be understood without go- 
ing through the derivation below. 

The derivation takes as starting point the expressions 
for the entanglement entropy 

S = -Tr{Mln A/ + (1 - M) ln(l - M)} (2.9) 

and the cumulant generating function 

Inx(A) = lndet{[l + (e'^ - l)M]e-^^Q} (2.10) 

= Trln{[l + (e^^ - l)M]e-'^'^} (2.11) 

in terms of the correlation matrix AI and a phase fac- 
tor e"*'^'^ which eventually drops out of the problem. 
For lattice models involving non-interacting fermions, the 
entries of M are the Green's functions Mij = (cijCii) of 
the fermionic degrees of freedom [64] as described in the 
Supplementary Material, while for a generic two-terminal 
scattering situation M is determined by the scattering 
matrix of the problem [65, 66]. Indeed, Eq. (2.11) is 
the famous Levitov-Lesovik determinant formula [47]. In 
both cases M has real eigenvalues in the range to 1 due 
to fermionic statistics. The similarity between Eq. (2.8) 
and Eq. (2.11) highlights the fundamental result that in a 
two-terminal (bipartite) setup charge transfer can be un- 
derstood in terms of simple elementary processes [65, 66]. 

To relate the von Neumann entropy to the cumulants 
of charge fluctuations, we first expand the logarithms of 
Eq. (2.9) around M = 0, 1 to obtain 



(2.12) 



with coefficients 



Ak = Tr[M(l - M)* + Af*^(l - M)]. (2.13) 

In the next step we relate Tr(Af'^) to the charge statis- 
tics described by Eq. (2.11). To this end the factorial 
cumulants (see Supplementary Material) 

Ffc = atlnx/(A)U=o (2.14) 

turn out to be particularly useful, where 

lnx/(A) =lnx(-iln(A + l)) (2.15) 

is the factorial cumulant generating function. The facto- 
rial cumulants corresponding to the generating function 
(2.11) are given by [50] 

i^fc = (-l)'-'(fc-l)![Tr(A/'=) 



Q] 



(2.16) 



Solving for Ty{M'') and substituting into Eq. (2.13) after 
having used the binomial theorem gives 



Ak = (-1) 



Fk , Fk 



(fe-1)! 



+ - 



fc! 



r=0 



k\ Fr+l 



(2.17) 



or, from Eq. (2.12), 



fc=i 



(-1)* 



k 



(fc-1)! 



"IT 



r=0 



rlk 



(2.18) 

The phase factor Q does not contribute to the sum. 

Eq. (2.18) expresses the von Neumann entropy purely 
in terms of factorial cumulants, which in principle are 
measurable quantities. It is nevertheless convenient, and 
somewhat more elegant, to express S directly in terms of 
ordinary cumulants. This is simply a matter of algebra, 
but we present here some important intermediate steps 
that would allow the interested reader to reproduce the 
derivation. 

The factorial cumulants can be expressed as a linear 
combination of the ordinary cumulants as 



Fk>i = Y,i-^)'~'"Siik,n)Cn, 



(2.19) 



where Si{k,n) are the unsigned Stirling numbers of the 
first kind. We can therefore write Eq. (2.17) as 



fe+i 

Ak = J2(-^)"~' 



k r+1 



r=0 n=l 



Si(k, n) Si{k + 1, n) 
(fc - 1)! fc! 

k\ Si{r + l,n) 



C„. (2.20) 



We note that Si(fc,fc + 1) = 0. Using the recursion re- 
lation [67] S'i(fc + 1, n) = <S'i(fc, n — 1) + kSi{k, n) in the 
first term and switching the order of sums in the second, 
we can simplify A^ to 



k+l 

Ak=Y.{{-l) 



Si{k,n — 1) 



fc! 



fc\ 5i(r + l,n) 



(2.21) 



Substituting into Eq. (2.12) gives (we can exclude n = 1 
since it contributes zero to the sum) 



k=l n=2 



S'i(fc, n — 1) 
fc!fc 

k\ Si{r + 1, n) 



rlk 



Cn 



(2.22) 



which we can write as a series in terms of the C„ by 
introducing a cutoff K in the outer sum and switching 
the order of the sums to get 



A'+l 



S= lim Va„(K)C„ 



(2.23) 
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where the cutoff-dependent coefficients an{K) are 

K 



k=n — l 
k 



EM 



Si{k,i 



1 



k\ 

k\ 5i(r + l,n) 
r 



(2.24) 



It can be sliown numerically for arbitrary integers n,k > 
n — 1 that the second term in square brackets is simply 



k\ Si{r + l,n) ^ ^ ^^„_i 5'i(fc,n- 1) 
rl r! A;! ' 



(2.25) 

so that we obtain Eq. (2.2). To the best of our knowledge 
Eq. (2.25) is a new identity for the Stirling numbers and 
a rigorous proof is desirable. As a check of the algebra, 
it is known that [68] 



Si(k,n- 1) 



k=n — l 



k\k 



(2.26) 



so that the coefficients (2.2) agree with the series pro- 
posed in Ref. [40] in the limit of infinite cutoff, K ^ oo. 
The similar series (2.3) for the Renyi is derived in Ap- 
pendix A. The convergence of both the von Neumann 
entanglement entropy (2.1) and Renyi entanglement en- 
tropies (2.3) is proved in Appendix B. 



B. Renyi Entropies and the Entanglement 
Spectrum 

As mentioned in the Introduction, many recent works 
[20, 21, 24, 25] have focused on the entanglement spec- 
trum, or the full set of eigenvalues of the reduced density 
matrix. The Renyi entropies of integer order are known 
to be fully equivalent to the entanglement spectrum [69] , 
which implies that, remarkably, for the non-interacting 
fermion systems considered here the charge statistics also 
encodes the entanglement spectrum. 

In this connection we present a very simple method for 
computing the full entanglement spectrum of a (finite) re- 
duced density matrix from the corresponding Renyi en- 
tropies of integer order. Specifically, for D degrees of 
freedom in the reduced density matrix only the Renyi 
entropies of order 2, . . . ,D are required. Although our 
interest here is on non-interacting fermions, the method 
is in fact completely general and is applicable to any den- 
sity matrix and hence of considerable interest to studies 
where the Renyi entropies but not the von Neumann en- 
tropy can be computed [43, 46]. The mathematical basis 
for this result, the Newton-Girard formulas [70], is long- 
established but to the best of our knowledge is not well- 
known in the physics literature. The method presented 
here should be compared with the standard method of 



analytically continuing the Renyi entropies iSa to a — > 1 
(see below) and the method presented in Ref. [69]. 

Suppose we are given & D x D reduced density matrix 
p. The Renyi entropies are defined for general a > as 



iSq — 



1 



1 - Q 



ln[Tr(p° 



(2.27) 



It is usually noted that the von Neumann entropy can be 
obtained by taking the limit as a — > 1, 



iS = — Tr(plnp) = lim Sa- 

a— >1 



(2.28) 



For convenience let 7?„ be the trace of the a-th power of 
the reduced density matrix. 



Rc = Tr(p' 



(2.29) 



Note that i?i = 1 for a properly normalized density ma- 
trix. 

Define the D x D matrix 



E ■ 



1 

R2 



Rd-1 
V Rd 



Rd- 
Rd- 



1 

R2 



\ 



D-1 
1 / 



(2.30) 



i.e., a quasi-lower triangular matrix with i?i = 1 on the 
main diagonal, R2 on the sub-diagonal, R,^ on the (ra— 1)- 
th sub-diagonal, 1, 2, 3, . . . , D — 1 on the super-diagonal, 
and zero everywhere else. For example, for D = 5 we 
have 



E ■ 



/I 1 0\ 
R2 I 2 00' 

i?3 i?2 1 3 

i?4 i?3 i?2 1 4 

\i?5 i?4 i?3 i?2 1/ 



(2.31) 



Denote the matrix obtained by taking the first n x n 
submatrix of E as En- Then we can form the polynomial 



P(x) = V^^(deti?Jx^- 



(2.32) 



n=0 



with the understanding that deti^o = 1- 

The polynomial P{x) is actually the characteristic 
polynomial of the density matrix p: P{x) = det(a;/ — p). 
If the individual eigenvalues of the reduced density ma- 
trix are desired then we can simply find all the roots of 
P{x). For example, in the simple case where the reduced 
density matrix represents a pure state the Renyi entropies 
would be zero and hence R2 = R3 = ■ ■ ■ = Rd = ^- Then 
(written here for the case D = 5) 



E : 



/I 1 0\ 
'11200' 
1113 
11114 
Vl 1 1 1 1/ 



(2.33) 
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FIG. 3: (color online). Spectral density function correspond- 
ing to tlie spectrum in Eq. (2.37), using Eq. (2.36) with 
e = 0.01. 

It can be shown that, aside from deti?o = Ij only 
det i5i = 1 is nonzero. The corresponding characteris- 
tic polynomial is therefore 

P{x)=x°-^{x-\), (2.34) 

so that the eigenvalues are 1 (multiplicity 1) and (mul- 
tiplicity D — 1) as expected. In contrast, for the fully 
mixed state where p has 1/D on the diagonal and zero 
everywhere else, we have i?„ = D^~" so that the charac- 
teristic polynomial is 

P{x)=(^x-^^ (2.35) 

and all of the eigenvalues are 1/D (multiplicity D) . The 
essential point, of course, is that we can obtain this not 
from knowing the full density matrix but from only the 
Renyi entropies and Eq. (2.32). 

For analytical purposes, it may be advantageous to 
write down the spectral density function 

^(z) = - lim Im lnP(z - ie), (2.36) 

where e is a small positive number that can be taken 
arbitrarily close to zero for increased accuracy. As an 
example, suppose the eigenvalues of the reduced density 
matrix p are 

0.0489694, 0.145703, 0.145703, 0.301666, 0.357959 

(2.37) 

where we have included a pair of degenerate eigenvalues 
to illustrate the effects of degeneracy. Then by perform- 
ing the above calculation for P{x) and using the resulting 
polynomial in Eq. (2.36) with e = 0.01, the spectral den- 
sity function appears as shown in Fig. 3. 

For a long time it has been possible to determine the 
von Neumann entanglement entropy only with exact di- 
agonalization (ED) for small systems and DMRG for one- 
dimensional and quasi-one-dimensional systems. For the 




FIG. 4; (color online). Quantum point contact (QPC) with 
chiral edge states (yellow lines) running along the edges of a 
sample containing a two-dimensional electron gas in a strong 
perpendicular magnetic field. The sample is divided into a 
left (L) and a right (R) region by a split gate acting as a 
QPC. Electrons in the two in-coming edge states are trans- 
mitted through the QPC with probability D or reflected with 
probability 1 — D. The number of electrons that have been 
transmitted, say, from the left to the right region, is denoted 
as m, which can be both positive and negative depending on 
the direction of the tunneling events. 



Renyi entanglement entropies of integer order n > 1 it 
was recently shown that they can also be determined in 
quantum Monte Carlo (QMC) calculations in two dimen- 
sions using an elegant trick [43, 46]. The von Neumann 
entanglement entropy, however, has remained inacces- 
sible outside of ED and DMRG, although it has been 
assumed that the Renyi entropies could be analytically 
continued to obtain the desired result. With the method 
presented here it is now possible to obtain not only the 
von Neumann entanglement entropy but also the full en- 
tanglement spectrum, assuming the Renyi entropies can 
be obtained to sufficiently high order. 



C. Applications 

We now illustrate the use of our formulas on a number 
of different examples. In one dimension, we expand on 
the analysis given in Ref. [36] of the entanglement en- 
tropy generated in a QPC under both equilibrium and 
non-equilibrium conditions, as well as of the spin-1/2 
XX model. For the latter, we discuss the influence of 
finite temperature and strong disorder in the chain. For 
non-interacting fermions in two dimensions we show how 
the cumulants reproduce the violation of the strict area 
laws for free fermions, while, in contrast, for bulk states 
of integer quantum Hall systems the area law is strictly 
obeyed. 
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1, Quantum Point Contact 

We first consider tlie entanglement entropy generated 
between two electronic reservoirs connected via a QPC, 
Fig. 4. In the situation discussed here, the chiral edge 
states of a two-dimensional electron gas in a strong per- 
pendicular magnetic field constitute the one-dimensional 
leads from which incoming electrons are either transmit- 
ted through the QPC and appear in the other reservoir 
or are reflected back into the out-going edge channel. We 
analyze the entanglement entropy between the two Fermi 
seas separated by the QPC. 

Without a bias difference between the electronic reser- 
voirs, the cumulant generating function for the number of 
charges m that have been transferred through the QPC 
(say, from left to right) is well-known and reads [71] 

Inx(A) = ^G(i) (2.38) 

with 

G(i) = In I ^sinh—1 (2.39) 

determining the time dependence. Here = h/{nkBT) 
is the time scale associated with the temperature T and 
Tc is an ultraviolet short-time cutoff. Additionally, the 
counting field A enters via the relation 

sin y = sin ^ (2.40) 

involving the QPC transmission D. It is assumed that 
the QPC is initially pinched off and instantaneously 
opened at time t = when the counting of transmit- 
ted electrons begins. A more detailed description of the 
opening of the QPC and its influence on the charge statis- 
tics can be addressed using the numerical procedure de- 
scribed in Ref. [72], but is not essential here. 

Interestingly, for perfect transmission D = 1, 
Eq. (2.40) immediately yields A = A*, implying that the 
fluctuations described by Eq. (2.38) are purely gaussian 
with zero mean, variance C2 = G/tt^, and all higher cu- 
mulants equal to zero. In this case, we can directly take 
the limit of infinite cutoff in Eq. (2.1) such that the coef- 
ficients an{K — > 00) = 2("(n) are given by the Riemann 
zeta function (cf. the generalization to the case with frac- 
tional quantum Hall reservoirs in Ref. [73]). Since only 
the second cumulant is nonzero, we find S = G/?> having 
used C(2) = 7r^/6. Considering now the zero tempera- 
ture limit — >■ 00 for which G(t) = ln(i/Tc), we readily 
obtain the simple result 

S = \\n-. (2.41) 

3 Tc 

Remarkably, this expression agrees with the prediction 

5=^ In- (2.42) 
3 a 



for a window in space of size £ of a conformal field the- 
ory with central charge c = 1 and short-distance cut- 
off Q [6, 7, 37, 74]. In our case, an analogous window 
[0, t] in time is used during which particles can delocal- 
ize among the two reservoirs, thereby making them en- 
tangled. It should be noted that this problem of local 
quantum quenches and the analogy between Eqs. (2.41) 
and (2.42) is qiute non-trivial and has been addressed 
extensively in Ref. [75]. 

In principle, the logarithmic dependence on time of the 
gaussian equilibrium quantum noise of a QPC is mea- 
surable. Indeed, measurement of the second charge cu- 
mulant showing a logarithmic dependence on time alone 
would constitute a strong indication of many-body en- 
tanglement. However, measurements of the equilibrium 
quantum noise should be contrasted with several other 
sources of noise that we discuss in turn below. Firstly, 
an imperfect transmission of the QPC, D < 1, leads to 
non-gaussian fluctuations which modify the result above. 
Additionally, any nonzero temperature causes classical 
thermal fluctuations, such that the von Neumann entropy 
eventually becomes linear in time. This is also the effect 
of a finite bias voltage, which gives rise to shot noise due 
to individual binomial charge transfer events that con- 
tribute only classically to the entropy increase without 
generating any entanglement between the reservoirs. 

In the case of imperfect transmission, the entanglement 
entropy maintains its logarithmic growth with time but 
the prefactor changes. This is illustrated in Fig. .5a, show- 
ing the entanglement entropy as a function of time for a 
QPC with transmission D = 0.5 obtained from the se- 
ries (2.1) with increasing cutoff number K. For K = 2, 
where only the second cumulant is included, we immedi- 
ately find S = [5/(27r)2]ln(t/rc) ~ 0.13 In(V'rc). This is 
already a good approximation of the exact result to which 
the series basically converges for K ~ 30. Importantly, 
the result obtained for K = 2 provides a lower bound on 
the full entanglement entropy, implying that a measure- 
ment of the second cumulant showing a logarithmic de- 
pendence on time would constitute a strong indication of 
many-body entanglement. For comparison, we also show 
Eq. (2.41) corresponding to a perfectly transmitting QPC 
with D = \. The entanglement entropy is maximal for 
D = 1, since the entanglement is generated by charge- 
neutral processes in which the two regions exchange par- 
ticles [40]. For imperfect QPCs where incoming electrons 
may refiect back on the QPC, the entanglement entropy 
decreases with decreasing D and, ultimately, when the 
QPC is completely pinched off with D = no electrons 
are transmitted and no entanglement is generated. 

We next discuss the inffuence of a nonzero electron 
temperature. Due to charge conservation, the even cu- 
mulants are the same if we count the number of electrons 
that have been transferred from right to left instead of 
from left to right. The von Neumann entropy of the 
two regions with respect to each other is therefore iden- 
tical even at finite temperatures, where the full system 
is no longer in a pure state. In Fig. 5b we show the 
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FIG. 5: (color online). Entanglement entropy in a quantum point contact (QPC). (a) Zero-temperature results for the time 
dependent entanglement entropy with different QPC transmissions D. The D = I result shown in red is given by Eq. (2.41). 
For D = 0.5, results were obtained from the series (2.1) with increasing cutoff number K from bottom to top. The thick blue 
line is the converged result for K = 30. The ultraviolet short-time cutoff is Tc = 10~^ ns. (b) Finite-temperature results for 
T = 10 mK (or r^a ~ 1.5 ns), D = 0.5, and Tc = 10~^ ns. Blue lines show results obtained with increasing cutoff number 
K and the thick blue line is the converged result. For comparison, zero and high temperature limits are indicated with a 
red and a blue dashed line, respectively. For short times t < Tfs the time-dependence is logarithmic, while for long time the 
high-temperature behavior eventually prevails and the entropy grows linearly with time, (c) Results for a biased QPC as a 
function of the transmission D. Here v = eV/{2kBT) is the ratio of the applied voltage V over temperature T. In the long-time 
limit t ^ T/}, the ratio 5(-D)/5(0.5) does not depend on time. The cutoff number is K = 10. As the voltage is increased, the 
entanglement entropy changes from a nearly linear dependence on D to that of a binomial event with success probability D, 
given by Eq. (2.45) and shown with a red dashed line. 



von Neumann entropy as a function of time at a realis- 
tic temperature of T = 10 mK and associated time scale 

— l.b ns. For short times t < Tp, the time dependence 
remains logarithmic since G{t) ~ ln(t/rc) and the be- 
havior is nearly identical to the zero temperature result 
indicated by a red line. In contrast, at long times t > Tp, 
temperature effects become important and the linear-in- 
time dependence G{t) ~ t/rp + ln[T^/(2Tc)] (shown as 
a dashed line) for high temperatures eventually prevails. 
Importantly, the zero-temperature behavior observed on 
sub-nanosecond time scales may be experimentally ac- 
cessible with current high-frequency noise measurement 
techniques that operate at bandwidths in the gigahertz 
regime [56, 58, 76-78]. 

Finally, it is instructive to consider the effect of a finite 
bias between the two Fermi seas. Our series are also 
applicable to such non-equilibrium situations. Focusing 
only on the thermal noise and the shot noise due to the 
applied bias, the cumulant generating function reads 

lnx(A) = -u+u_, (2.43) 

where 

u± = V ± cosh""'^ {D cosh(u + i\) + {\ — D) cosh(7;)} 

(2.44) 

and V = eV/{2kBT) is the ratio of the applied voltage 
over temperature [71]. Fig. 5c shows the von Neumann 
entropy as a function of the QPC transmission D for 
several applied voltage differences. For low voltages, the 
entropy generated is mainly due to thermal fluctuations 
and depends nearly linearly on D. As the bias voltage 



is increased, shot noise-dominated processes become in- 
creasingly important and in the regime, where the ap- 
plied bias voltage is much larger than temperature, the 
generated entropy is solely due to individual binomial 
charge transfer events occurring with probability D. In 
that case, the von Neumann entropy is that of a classical 
binomial process involving eVt/h attempts for which 

eVt 

S = —[DlnD+{l-D)\n{l-D)]. (2.45) 

This expression is shown with a dashed red line and is 
in excellent agreement with our results for the high-bias 
limit. Eq. (2.45) is also the electron-hole entanglement 
predicted for a biased tunnel junction [79], and in both 
cases the Renyi entropies for small D are 

T) 

Sn = tN, (2.46) 

n — I 

where A'^ = {eVt/h)D is the number of transferred parti- 
cles [80] . We note that in the finite- bias situation consid- 
ered in the last part, the generated von Neumann entropy 
is only due to classical fluctuations and does not reflect 
true many-body entanglement. 

2. The Spin- 1/2 XX Chain 

In our next example, we consider the one-dimensional 
spin- 1/2 XX model. We first analyze the situation where 
the coupling between neighboring spins is constant along 
the chain and discuss in this case the influence of a finite 
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temperature on the spin fluctuations and the entangle- 
ment entropy in the chain. Secondly, we turn to the situ- 
ation where the coupling between neighboring spins has 
random magnitude. Importantly, the relation (2.1) be- 
tween the von Neumann entanglement entropy and the 
cumulants is linear, such that the equation also holds 
when averaged over the disorder. The series then be- 
comes a linear relation between the disorder-averaged 
entanglement entropy and the disorder-averaged cumu- 
lants. 

The spin-1/2 XX model in one dimension is described 
by the Hamiltonian 



(2.47) 



where Sf, a = x,y,z, is the operator for the a- 
component of the spin on site i. We first consider Ji = J 
being constant along the chain and discuss analytically 
the zero-temperature and high-temperature behaviors of 
both the von Neumann entanglement entropy and the cu- 
mulants of the spin fluctuations. These results serve as 
important limiting check points for our numerical calcu- 
lation with intermediate temperatures. Unless otherwise 
stated, in this section all system sizes count the number 
of lattice sites, i.e., lengths are in terms of the lattice 
spacing and are dimensionless. 

The zero-temperature result for the entanglement en- 
tropy was determined by Calabrese and Cardy who found 
[7] 



5(L/2,L;r = 0) 



1, L 

-In - + si, 

3 TT 



(2.48) 



where L is the total length of the spin chain (in units 
of the lattice spacing) and the subsystem is taken to be 
half the chain of length L/2. The constant si was ob- 
tained analytically by Jin and Korepin and can be found 
in Ref. [81]. 

Although it is difficult to compute the full set of cumu- 
lants analytically, it is worth obtaining analytical results 
for the fluctuations for the spin-1/2 XX chain [53]. For 
periodic boundary conditions (PBCs) the computation 
in the Supplementary Material gives 



(2.49) 

plus 0{£~'^) corrections, with /i = 1 + 7 + In 2. For 
open boundary conditions (OBCs) the computation in 
the Supplementary Material leads to 

-^x°x-(^)4-F-C(2.) + ^ 

ln(2^) + 7 + ln2] 
■[ln(2i') - ln2] (2.50) 



7r2(2^) 
7r2(2^): 



plus 0{£~^) corrections. In both cases the linear term 
vanishes, and for the case of OBCs subleading oscillations 
exist. 



It is interesting to derive Eq. (2.49) for PBCs in a way 
difl'erent from the computation presented in the Supple- 
mentary Material that illustrates the similarity to the 
computation of the entanglement entropy for the spin- 
1/2 XX chain in Ref. [81], although it is less controlled 
than the previous method. The method also suggests 
how one might be able to compute all of the cumulants. 
In principle, the full set of cumulants can be calculated 
from the generating function 



x(A) = (exp(jA^S'| 



exp I i\ fijaj 



(2.51) 
(2.52) 



where a|, cij are the Jordan- Wigner fermion creation and 
annihilation operators (cf. Supplementary Material) . We 
have 



X(A) = det 



cos^j/+ (jsin^ |M 



(2.53) 



where M is the £ x C matrix with elements 



2(_l^(^-J + l)/2_J_ fori-jodd. 







for i — j even. 



(2.54) 



For f — > CO the determinant (2.53) can be evaluated using 
the generalized Fisher-Hartwig conjecture for Toeplitz 
determinants to obtain the result (we assume A < tt since 
we desire the derivatives at A = for the cumulants) [82] 



liH(A) = -A_i„(2t) + 21n 



where G(z) is the Barnes-G function defined by 
Gil + z) = (27r)^/2g^p 



27r 
(2.55) 



1 1/, ^ 2 
-z (1 + -f)z^ 

2 2^ " 



n 



1+ 



exp 



z + 



2k 



(2.56) 



Using 

ln[G(l + z)G{l - ^)] = -(1 + 7)2^ + 0[z^) (2.57) 

we find that ^xx^ = C2 = (— ^^a)^ lnx(A)|A=o agrees 
with Eq. (2. 19). 

At finite temperature the von Neumann entropy can 
be found from the correlation matrix M of the full chain 
given in the Supplementary Material (see also Ref. [83]. 
In the high-temperature limit /3J <C 1 the correlation 
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matrix is 



^=2 



1 

-/3J/4 



\-p/?J/4 



-/3J/4 
1 



-/3J/4 



1 

-/?J/4 



-/3J/4 
1 



(2.58) 

with p determined by the choice of boundary conditions 
(Supplementary Material); here we restrict to p = ±1 
for PBCs. Considering a subsystem of length £ < L, 
the eigenvalues of the correlation matrix restricted to the 
subsystem A/'^' are easily found to be 



1 



1 



(3 J Tin 

— cos 

2 £+ 1 



(2.59) 



for n = !,...,£. Using the expansion of the binary en- 
tropy function H2{x) = —x\nx — {\ — x) ln(l — x) 

H,(^-^^=\n2-'^+0{e^) (2.60) 

for e <C 1 and the identity 



the von Neumann entropy (2.9) becomes 



5(AL;/3J«; 1) 



4 



ln2 



PI 
4 



(2.61) 



I. (2.62) 



For comparison, if we consider the full chain i = Lwe can 
similarly demonstrate that the von Neumann entropy is 
given by 



5(L,L;/3J<C1) 



ln2- 



/3J 
4 



(2.63) 



Eq. (2.63) is of course nothing but the thermodynamic 
entropy of the XX model at high temperature, which is 
purely extensive in the system size L. Eq. (2.62), mean- 
while, is almost extensive in £, and the constant contri- 
bution becomes negligible with increasing temperature, 
i.e., as PJ decreases. This is in contrast to the zero- 
temperature behavior (2.48) where the logarithmic de- 
pendence on the subsystem size is a clear signature of 
true quantum correlations, i.e., entanglement. 

We now compute the cumulants C„ for the total spin 
5*^ in the subsystem. Since we know the eigenvalues of 
M^^' we can write for £ < L 



Inx(A) = ^ ln{[l + (e^^ - l)i.„]e-V2} 



(2.64) 



In cos 



A 
2 



i tan — I cos 



PJ 
2 



Trn 
£ + 1 



+ - tan^ - 



(3 J nn 

— cos 

2 £+1 



+ 0{{l3jf). 
(2.65) 



The second term of Eq. (2.65) is an odd function of A 
and does not contribute to the even cumulants, so we 
may write the effective cumulant generating function to 
leading order in (3 J as 



Inx(A) 



I3J 
4 



tan 



A 



1 ^^fPJ 
m cos — + 



tan 



Similarly, for £ = L 
Inx(A) 



In cos — + 



tan 



L. 



(2.66) 



(2.67) 



It is easily shown numerically that the cumulants aris- 
ing from Eqs. (2.66), (2.67) lead to the entanglement en- 
tropies given in Eqs. (2.62), (2.63), respectively. 

In Fig. 6 we illustrate how the series (2.1) converges to 
the von Neumann entropy at various temperatures as the 
number of cumulants is increased. Exact results were ob- 
tained by numerically diagonalizing the correlation ma- 
trix (see Ref. [84] and the Supplementary Material). The 
results interpolate nicely between the zero and high tem- 
perature limits given by Eqs. (2.48) and (2.62), respec- 
tively. In particular, at low temperature where the en- 
tanglement entropy is similar to that of a pure state, the 
entanglement is essentially equal for the subsystem and 
its complement, such that S{£, L) ~ S{L — £, L) as seen 
in Fig. 6a. In the high-temperature regime. Fig. 6d, this 
symmetry is broken and the von Neumann entropy for 
the subsystem and its complement are now very differ- 
ent. This is clearly an effect of the finite temperature, 
which causes random and independent spin flips on each 
site. This effect is extensive in the subsystem size £ and 
causes the even cumulants of the spin fluctuations in the 
subsystem and its compliment to be different. 

For subsystem size £ = L/2 the temperature T* where 
thermal effects become important can be estimated from 
CFT arguments as follows. For an infinite system and 
subsystem size x at finite temperature T = 1/(3, we can 
use the mapping [7] 



Vl3 TTX 

— sinh — 

TT Vj3 



(2.68) 



where v is the effective velocity, in this case v = J. The 
sinh becomes important when nx ~ vf3, or with x = L/2 
[83], 



2J_ 



(2.69) 



This can already be seen in Figs. 6b and 6c with temper- 
atures around T* . For T = 2T* , Fig. 6c, the von Neu- 
mann entropy is clearly not symmetric. The crossover 
between zero-temperature and high-temperature behav- 
ior is further illustrated in Fig. 7 showing the entangle- 
ment entropy as a function of temperature for differ- 
ent system sizes L. Here the partition is taken in the 
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FIG. 6: (color online). Exact entanglement entropy (dashed) 
and approximation by cumulants (solid lines) of the spin-1/2 
XX chain, L = 100, as a function of subsystem size, at temper- 
ature (a) T = 0.5r*, (b) T = T*, (c) T = 2T\ (d) T = 50T* 
where T* = 2J/{itL) is the temperature at which thermal ef- 
fects become significant. The cutoff K increases from bottom 
to top. Note the differences in scale on the vertical axes. 
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FIG. 7: (color online). Entanglement entropy of the spin-1/2 
XX chain as a function of temperature T, for several system 
sizes L. The partition is taken in the middle at £ = L/2. 
Solid lines show analytical results for the zero- and high- 
temperature limits, and dashed lines indicate the crossover 
temperature T* , Eq. (2.69), where the entanglement entropy 
begins to deviate significantly from the zero-temperature be- 
havior. 



middle of the chain and periodic boundary conditions 
have been imposed. The solid lines show the analyti- 
cal zero- and high-temperature limits, while the vertical 
dashed lines correspond to the crossover temperature T* 
defined by Eq. (2.69). The crossover temperature T* sep- 
arates the von Neumann entropy into a low-temperature 
regime, governed by zero-temperature behavior with log- 
arithmic dependence on the subsystem size i, and a high- 
temperature limit, where the von Neumann entropy de- 
pends linearly on i. In the low-temperature regime, the 
generated von Neumann entropy is mostly due to quan- 
tum correlations or entanglement, whereas in the high- 
temperature regime the von Neumann entropy mainly 
reflects classical, thermal fluctuations. This is important 
for actual experiments, where one wishes to measure the 
entanglement and not classical fluctuations. 

Finally, we briefly discuss the effect of disorder. We 
thus allow the coupling strength Ji between sites i and 
i + 1 to vary along the chain with probability distribu- 
tion P{Ji). For essentially any distribution P{Ji), the 
fixed point of the model is known to be described by the 
random singlet phase (RSP) in which the ground state 
is a product of singlets between spins at arbitrarily large 
distances [85]. Each singlet = (|t4j) - |iitj))/v^ 
that crosses the boundary between a subsystem of length 
I and the remainder of the system contributes In 2 to the 
entanglement entropy. It is thus clear that we can write 
[74] 



>{i) = filn2. 



(2.70) 



where n is the disorder-averaged number of singlets that 
cross the boundary between the two blocks. It has been 
shown that n ~ (1/3) In^ [74, 86, 87]. 
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For a factorized state like the RSP, it is also easy to 
find the cumulants of the total spin in the subsystem. 
For a subsystem of length £ the generating function for 
S*^ is 



x(A) = (*RSp|exp(iA^^| 



I* 



RSP/ 



(2.71) 



where |3'rsp) = Y\.{ij)\^j) ^ product of non- 
overlapping singlets covering the entire system. Singlets 
between two spins within the subsystem ( "in- in" ) or two 
spins outside the subsystem ("out-out") are eigenstates 
of Sf + S'j with eigenvalue and contribute a factor of 
e° = 1 to the expectation value. In contrast, a singlet be- 
tween one spin inside the subsystem and one outside ( "in- 
out") contributes exp[iA(S'f + S^)] \ij) = cos(A/2) to 
the product, so if there are n in-out singlets then 



X(A) 



(2.72) 



The disorder-averaged cumulant generating function is 
therefore 



Inx(A) = n In cos 



A 



(2.73) 



Eq. (2.1) applied to Eq. (2.73) then leads to Eq. (2.70). 
As can be shown using exact diagonalization, moreover, 
the formula holds exactly for each separate realization of 
the disorder so that the equivalence between the cumu- 
lants and the entanglement entropy is a true equivalence. 

It is worth noting that the previous discussion does 
not apply if the disorder is locally correlated [88]. In this 
case the fixed point is not a RSP and the prefactor of 
the logarithmic growth of the entanglement entropy S 
with subsystem size is not universal. Since the model is 
still that of non-interacting fermions, however, Eq. (2.1) 
still holds. Perhaps more suggestively [cf. Eq. (3.33)], 
S/C2 ~ 7r^/3 for locally correlated disorder. 



3. Entanglement Entropy of Fermions in Two Dimensions 

We now consider systems in dimension higher than 
one and focus in particular on the area laws obeyed by 
systems in two dimensions. It has already been well- 
established through explicit arguments [89, 90] and CFT 
arguments on the Fermi surface [12] that for free fermions 
the entanglement entropy of a subsystem of linear size L 
in d-dimensions typically obeys an area law with multi- 
plicative logarithmic corrections: 



S ~ L'^-MnL, 



quasi-area law. 



(2.74) 



For the particular case of free fermions on a square 
lattice with PBCs at half-filling as in Fig. 8, a numerical 
fit shows that the prefactor of the leading L In L term is 
almost certainly exactly 1/3 as predicted by a formula 




FIG. 8: (color online). Entanglement entropy and approxi- 
mation by cumulants for free fermions on an L x L square 
lattice with PBCs, where the partition is taken to be half 
the system. The solid line shows a fit to the generic form 
S = aLlnL + bL + c\nL + d + e/L, where, in particular, 
a ~ 0.333. 



based on the Widom conjecture [90, 91]. Of course, un- 
like the one-dimensional case the prefactor is not univer- 
sal and depends, for example, on the chemical potential. 
The behavior of the cumulants is necessarily similar as 
illustrated in Fig. 8. The approximation of the exact en- 
tanglement entropy using the series (2.1) improves, and 
we eventually find the quasi-area law (2.74) for the entan- 
glement entropy with prefactor 0.333 in agreement with 
Ref. [91]. The similar behavior of the cumulants can be 
understood from the CFT nature of the problem as we 
discuss in further detail in Sec. Ill A 4. 

Of course, there also exist systems of non-interacting 
fermions in d dimensions which obey the usual strict area 
law 



d-l 



strict area law, 



(2.75) 



without multiplicative logarithmic corrections. Ref. [91] 
has explored this behavior as a function of the dimen- 
sionality and shape of the Fermi surface, showing that 
entanglement entropy (and hence the charge statistics) 
contains important information about the geometry of 
the Fermi surface. The argument of Ref. [12] elegantly 
illustrates the role of the Fermi surface in the appearance 
of the multiplicative logarithmic correction for Fermi 
liquids: one essentially breaks the Fermi surface into 
0{L'^~^) "patches," with each patch described by a c = 1 
CFT and hence contributing (1/3) In L to the entangle- 
ment entropy. 

Here we consider an important example of a two- 
dimensional system of fermions that obeys a strict area 
law: the integer quantum Hall effect (IQHE) at filling fac- 
tor t/ = 1. The real-space entanglement entropy for this 
system was computed for various geometries in Ref. [92] 
and found to obey a strict area law. Following Ref. [92] , 
we thus consider the lowest Landau level of electrons in a 
cylinder of size x Ly in the gauge A = 5(0, x) where 
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FIG. 9: (color online). Cylinder geometry for the integer 
quantum Hall effect. The "area" of subsystem A is 2Ly. 



B is the magnetic field, i.e., the magnetic field is per- 
pendicular to the surface of the cylinder. The cylinder is 
periodic in the y-direction and infinite in the x-direction. 
We will work in units where the magnetic length is equal 
to 1 [92]. In the limit L^^Ly ^ 1 the ground-state cor- 
relation matrix Mrr' = {c\ct'), where Cr is the fermion 
annihilation operator at position r, is given by 



Mr- 



1 

27r 



exp - ^ {x-x'f-^ {y-y'f-^ {x+x') iy-y') 

(2.76) 

We are interested in the entanglement entropy of the re- 
gion A defined by 



-f<a^<f, 0<y<Ly (2.77) 

as shown in Fig. 9. It was found in Ref. [92] that for 
2> 1 the eigenvalues of M^r' are 



Erfj^+M-Erfj'/i-l 



(2.78) 



where /i is a continuous index ranging from — oo to oo 
with measure Ly/{2Ti) and Erf(x) = (2/v^) dt e~* 
is the error function. Each eigenvalue contributes 
Lj,_ff2(r;(/i, i'a;))/(27r) to the entanglement entropy, where 
H2{x) = —X \nx — {1 — x) ln(l — x) is the binary entropy 
function. The total entanglement entropy is then given 
by the area law 



S{e,,Ly) = 2a{i,)Ly, 



with 



a (4) 



1 

An 



(2.80) 



In this geometry the "area" is given by the perimeter 
length 2Ly , so the result is expected to depend only on Ly 
for long cylinders 4^1- Indeed, in this limit the pref- 
actor a{£x) approaches the constant value ao — 0.203. 



FIG. 10: (color online). The exact (dashed line) entanglement 
entropy, under the assumptions that lead to Eq. (2.79), and 
approximation by cumulants for different cutoffs K, increas- 
ing from bottom to top, for the integer quantum Hall effect 
on a cylinder at filling factor u = 1. 



The factor of 2 in the area law (2.79) is due to the par- 
ticular choice of geometry, where the subsystem has two 
boundaries of length Ly aX x = ±1^/2 connecting it to 
its complement. For a subsystem consisting of a semi- 
infinite cylinder with only a single boundary of length Ly , 
the factor of 2 in Eq. (2.79) would be replaced by unity. 
The absence of a constant term, moreover, indicates the 
absence of topological order [92]. The "strip" geometry 
has also recently been used to study entanglement and 
Renyi entropies in the spin- 1/2 Heisenberg model on a 
square lattice [46]. 

These results also follow directly from the fiuctuations 
of charge across the boundaries between the subsystem 
and its compliment. To illustrate this, we evaluate the 
generating function for charge fluctuations (2.11) using 
the correlation matrix in Eq. (2.76). We find 



lnx(A,4) 



L 



d^i ln[l + (e^^- 1)77(^,4)], (2.81) 



with corresponding cumulants 



(7„(4) = 2/3„(4)i, 



(2.82) 



where 



(2.79) /3„(4) 



1 

47r 



d/i {-idxT ln[l + (e^^ - IHp, 4)] I A=o. 

(2.83) 

Fig. 10 shows the approximation by cumulants to the 
entanglement entropy [specifically, the prefactor q(£x) 
in Eq. (2.80)], plotted as a function of ^4 following 
Ref. [92] . As mentioned previously, the entanglement en- 
tropy (and also the cumulants) are independent of 4 for 
large i^, and depends only on Ly. We note that similar 
calculations can be carried out for the disk and spherical 
geometries [92]. 



15 



III. FLUCTUATIONS IN ONE-DIMENSIONAL 
INTERACTING SYSTEMS 



wavelength (|fc| <g; po, Po is the mean density) part of the 
density fluctuations in a LL is given by 



As shown in Sec. II, for systems that can be mapped 
to a problem of non-interacting fermions the complete 
set of even cumulants of the charge fluctuations fully 
determines the von Neumann and Renyi entanglement 
entropies, as well as the entanglement spectrum. This 
equivalence appears to be unique to non-interacting 
fermions [73] . For many one-dimensional conformally in- 
variant critical systems with a conserved U(l) charge, 
however, the fluctuations can be shown to scale logarith- 
mically with the system size, just as for entanglement 
entropy. We can also determine the prefactor which de- 
pends not on the central charge of the theory but on the 
compressibility (susceptibility) of the U(l) charge. In the 
following sections we show how this comes about and dis- 
cuss when and how the behavior of the fluctuations and 
entanglement entropy differ. The area-law behavior of 
fluctuations in gapped systems is also discussed. 

In the remainder of the paper we will be concerned 
exclusively with the second cumulant C2, which we refer 
to as the "fluctuations" T . 



A. One-Dimensional Critical Models 

The low-energy properties of many one-dimensional 
systems are captured by Luttinger liquid (LL) theory 
[93-95] and are described by the Hamiltonian 



LL 



V 

2tt 



dx 



(3.1) 



where v and K are the effective velocity and Luttinger 
parameter, respectively. We recall that for such a LL, a 
free bosonic CFT with central charge c = 1, the entan- 
glement entropy of a region of size x in an inflnite system 
scales as 



5LL(^) = ^ln^+Sl 



(3.2) 



where si is a non-universal constant and is sometimes 
combined into the logarithm term through a short- 
distance cutoff a as ln(^/a). Eq. (3.2) is for a system with 
PBCs; for OBCs the prefactor is divided by 2, 1/3 1/6, 
while (. — > 2£. For systems with flnite size L, moreover, 
conformal mapping leads to the simple substitution 



i- 



L 7r£ 
— sin — . 

TT L 



(3.3) 



1. Luttinger Liquids 

We now consider the fluctuations. In the following we 
measure distances in physical units and thus restore the 
short-distance cutoff to emphasize the ultraviolet diver- 
gence of the fluctuations and its regularization. The long 



p{x) 



^ Po + -9x0(.x) 

TT 



(3.4) 



so that for a block of length I extending from x = to 
X = £ we have 

1 , 



Na-{Na) = -Wx)-m], 



(3.5) 



where we have used the fact that (4>{x) — 4>{y)) = 0. The 
field (f) has the mode expansion 



0(a;) = 00 + 



n\g\/2 




(e^«"&, e-'^^&t) (3.6) 



where a is a short-distance cutoff, not necessarily the 
same as a for the entanglement entropy, and are 
bosonic creation and annihilation operators. For L — > cx) 
we can convert the momentum sum into an integral to 
find 



{m-mf) = 2K dq 



K 



In 



Q 
q2 



Thus for ^ > a 



7r2j-LL(£) = {[m - mf 



Kin-, 
a 



(3.7) 
(3.8) 

(3.9) 



As for the entanglement entropy, the same result with 
K — > K/2 and ^ — >• 2£ is obtained when there is a bound- 
ary, due to the constraint (j){0) = constant. This can be 
seen from 



0(2:) = 00 + « 

for OBCs. Then 

([0(£)-0(O)]2) 



9>0 



-aq/2 



sin(g.x)(&, - St) (3.10) 



-ctq 



K 

Iq ' q 



= — In 
4 



dq sin (g£) 

(2^)2 + a2 



Assuming £ ^ a we obtain 



.2J-OBC, 



— In — . 

2 a 



(3.11) 
(3.12) 

(3.13) 



We see that for both PBCs and OBCs in the limit 
^ — > cx) the leading terms of the entanglement entropy 
and fluctuations obey the relation 

S{€) ^ /^^.F(£), (3.14) 



which generalizes the prefactor 2(^(2) = 7r2/3 for free 
fermions [36, 40]. As a simple example, for the fermionic 
Hubbard model with spin-charge separation into two LLs 
[94] the ratio of the entanglement entropy to the fluctua- 
tions is controlled by for the charge degree of freedom 
and Kg = \ for the spin degree of freedom. 
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2. Case of Disjoint Intervals 

For the entanglement entropy, there has been growing 
interest in the case where the subsystem A is composed 
of disjoint intervals [7, 96, 97]. In particular, consider the 
case where A extends from xi to X2 and also to 2:4 (we 
consider PBCs once more). It was originally predicted in 
Ref. [7] that for a free boson the entanglement entropy is 
given by 

X ? 1 , Xi2X^iXi4,X22, , „ 

S[xi,X2,X3,Xi) = -in |-2si, (3.15) 

3 Xio,X2A 

cf. Eq. (3.2). It has been shown [96, 97] that this is not 
quite correct, due to the fact that for disjoint intervals the 
topology of the Riemann surface used in the derivation is 
highly non-trivial, i.e., non-local. It is interesting that, as 
we show below, the particle number fluctuations within 
LL theory, i.e., for a free boson, retains the simpler form 
of Eq. (3.15). This is a case where the entanglement en- 
tropy is clearly more sensitive to the full operator content 
of a GET. 

Eor the two-interval case considered above the number 
operator for subsystem A is 

Na-{Na) = -{[^(X2)-0(X1)] + [,^(X4)-<^(X3)]}, (3.16) 

SO that 

'k'^T-LI,{xi,X2,Xz,Xa) 

~ iv In 7, — 



+ 2^X2) - <f>ix,)MxA) - 0(X3)]) (3.17) 



= Kin 



X12X34 



K 



dq 



-a\q\ 



A' In 



■^123:343:143^23 
3:133:240^ 



(3.18) 

(3.19) 



where Xij = Xj — Xi. Thus Eq. (3.19) is essentially iden- 
tical to the naive result for the entanglement entropy in 
Eq. (3.15). 



3. Subleading Corrections 



For — 1 < A < 1 Hamiltonian (3.20) is gapless and de- 
scribed by LL theory. The isotropic point A = 1 is some- 
what special due to the effects of a marginally relevant 
operator and its infrared behavior is best considered in 
the context of the Haldane-Shastry model described be- 
low, while A = — 1 is not conformally invariant because 
of the quadratic spectrum. At A = the model reduces 
to the exactly solvable spin- 1/2 XX chain where the exact 
results of Sec. 11 apply, but it is instructive to analytically 
compute Sec. II C 2. 

From LL theory, for x ^ a the spin-spin correlation 
function in the ground state of Hamiltonian (3.20) is 
given by 



TT 

r,2m.^K-2 



Am cos(7rmr) 



771 > 2 



r 



27772 A' 



(3.21) 



where from the Bethe ansatz solution of the XXZ model 

-1 



1/ cos:^ 

2 V TT 



(3.22) 



while the Am are non-universal coefficients (interest- 
ingly, there also exists a conjectured formula for Ai, see 
Refs. [98, 99]). For the XX model, for example, K = I 
and Ai = 1/(271^), Am>2 = 0. Strictly speaking, of 
course, Eq. (3.21) is an asymptotic expansion. However, 
we will see that summing this expression gives results in 
excellent agreement with the numerical results, neglect- 
ing the linear term arising from the short-distance physics 
not taken into account within LL theory. With precise 
knowledge of the microscopic physics this linear term 
vanishes, as illustrated by calculations for the Haldane- 
Shastry chains described below and detailed in the Sup- 
plementary Material. It is also important to note that 
Eq. (3.21) is not correct for the isotropic point A = 1 
where marginal operators induce multiplicative logarith- 
mic corrections to correlation functions. 

To find the fluctuations, we first rewrite Eq. (1.3) as 



ijeA 



(3.23) 



As shown in Ref. [53] and the Supplementary Material, 
the fluctuations corresponding to Eq. (3.21) are then 
given by 



It is possible to compute subleading corrections to the 
fluctuations within LL theory. To illustrate, we consider 
the spin-1/2 XXZ model with PBCs at zero magnetic 
fleld (i.e., (S'f ) = 0), which is described by the Hamilto- 
nian 

Hxxz = J ^(Sf Sf+i + Sf Sf+i + A^f Sf+i)- (3.20) 



7r^Txxz{£) = Kln£ + f2-A, 



(-1)^ 

g2K 



(3.24) 



plus 0{£~^) corrections. This provides an alternative 
method for flnding the LL parameter K in numerical cal- 
culations instead of using the spin-spin correlation func- 
tion [53]. We note that for ii' < 1 the r"^^ term dom- 
inates at long distances over the Fermi liquid-like term 



17 



proportional to r~^. The logarithmic divergence, in con- 
trast, originates from the universal Fermi liquid- like term 
which decays as r~^, i.e., from the short-distance corre- 
lations. This is of course an important feature of the 
entanglement entropy, whose divergence is due to short- 
distance correlations at the boundary between two sub- 
systems. 

Since the SU(2)-symmetric point A = 1 is an impor- 
tant example of a model where the full symmetry group 
is actually larger than, but contains, U(l), it is worth 
checking that the above considerations are also valid for 
a closely-related system known as the Haldane-Shastry 
model [100, 101]. This model is in the same universality 
class as the Heisenberg model and has the peculiar fea- 
ture that its finite-size behavior is the same as its ther- 
modynamic limit, so that all logarithmic corrections are 
absent. The Haldane-Shastry model is described by the 
Hamiltonian 



(3.25) 



where Jy = l/d{i — j)^ and d{i,j) is the chordal distance 

TTX 



d{x) = - 
7r 



L 



(3.26) 



The exact spin-spin correlation is known to be [102] 

{sUrSt) - {sUr){S!) = h-iy^^, (3.27) 

where 



Si(x) = / dt 



sini 



(3.28) 



We note that the asymptotic expansion of Eq. (3.27) 
yields the same two leading terms as the LL expression 
in Eq. (3.21) with K = 1/2, but the exponents differ for 
the remaining terms. As shown in the Supplementary 
Material the fluctuations are given by 

^'-FHS = ;Tln£ + /HS-rT7^^ (3.29) 
2 lb I 

plus 0{E~'^) corrections, where /hs is related to an inte- 
gral whose numerical value is /hs/""^ — 0.197217. Fig. 11 
shows the excellent agreement between DMRG results 
and the analytical expressions, with the usual mapping 
£ — ^ (L/tt) sin(7r^/L) for finite systems. The calculation 
in the Supplementary Material also illustrates how ex- 
plicitly demonstrating the vanishing of the linear term is 
non-trivial and depends very much on the exact micro- 
scopic details of the theory. Since the vanishing of the 
linear, "volume" term is nevertheless universal, it is in- 
teresting that the short-distance behavior is not arbitrary 
but rather constrained by this necessity. The example of 
the Haldane-Shastry chain also shows that the leading 
logarithmic scaling of the fiuctuations can persist even 
when long-range interactions are introduced. 



0.35- 




0.25 



100 



FIG. 11: (color online). Fluctuations of the Haldane-Shastry 
model, L = 100. Numerical results from DMRG (squares) are 
in excellent agreement with the analytical result (solid line) 
in Eq. (3.29) with finite-size mapping. 



4- Extension to CFTs with a Conserved U(l) Charge 

For critical one-dimensional systems described by 
CFT, a more general argument may be given for the 
scaling of the fluctuations that does not rely on tradi- 
tional LL theory. Briefly reviewing the argument given 
in Ref. [53], this is because a conserved U(l) charge is al- 
ways described by a free bosonic generator in the theory, 
and consequently the generating function is given by 



X(A) 



(3.30) 



Then tt'^Ta = ii'^{-idxf\nxW\\=a = ffln(a;/a). The 
prefactor g can always be fixed by considering the phys- 
ical meaning of the conserved charge, but we can give a 
heuristic argument for its value as follows. At finite tem- 
perature 1//3 (we will set the effective velocity v = \ for 
now) the mapping z — > = [/3/(27r)] In z in Eq. (3.30) 
gives 



TT J-{x, p) = g In ( sinh 



(3.31) 



For sufficiently large x 3> /3 such that interactions across 
the boundary can be neglected (which is possible since 
correlations decay exponentially), we may consider the 
subsystem ^ to be a grand canonical ensemble in equi- 
librium with a bath consisting of the remainder of the 
system [103] . This is of course only possible if the total 
particle number is fixed. Then from standard statistical 
mechanics one has J-{x,l3) ~ Kx/fi where k = dn/dfi is 
the compressibility (susceptibility x = dm/dB for spins), 
so that by matching Eq. (3.31) for x ^ (],a we find 



g ■ 



(3.32) 



We have put in the velocity v for completeness. Note that 
this is consistent with the LL expression since K = ttvk. 



18 



Thus for such systems we observe that 



S{(.)r.—\T{t), 

nvK 6 



(3.33) 



which again generahzes the free-fermionic prefactor 
2C(2) = 71-2/3 [36, 40]. 

An important class of models is described by the 
SU(2)-symmetric Wess-Zumino-Witten non-linear a- 
model with topological coupling k [SU(2)fc WZW model]. 
They are CFTs with central charge c = 3k/{k + 2), so 
that the SU(2)i WZW model, with c = 1, is simply a 
free-boson theory. Furthermore, the spin susceptibility 
is given by nvx = k/2 where v is the effective velocity 
[53, 104]. A microscopic reahzation of the SU(2)fc WZW 
model for general k was found by Takhtajan [105] and 
Babujian [106] (TB) governed by the Hamiltonian 



(") 



(3.34) 



where P^'^^ is a projector onto the subspace with total 
spin V oi a. bond (i, i -\- 1). The TB Hamiltonian is solu- 
ble by Bethe ansatz and k = 25. We performed DMRG 
calculations of the Hamiltonian (3.34) for k = 2 and 3 
(fc = 1 being the spin- 1/2 Heisenberg chain) and found 
the predicted result nvx = k/2 within a considerable er- 
ror since the susceptibility is renormalized (like for the 
spin-1/2 Heisenberg chain). For S > 1/2 there are no 
models analogous to the Haldane-Shastry chain which do 
not have logarithmic corrections and renormalization ef- 
fects (see below) . The central charges are not affected by 
the renormalization and we obtained the expected values 
c = 3/2 and c = 9/5, respectively. Since k = 2S for the 
TB models, one might ask how we can in general know 
that we really measure k/2 and not just the spin S. In 
fact we can prove that we measure k/2: the half-integer 
spin-S' Heisenberg models with nearest-neighbor interac- 
tions are all known to be described by the same SU(2)i 
WZW model. That is, from the CFT point of view the 
5* = 1/2, 3/2, 5/2, etc. Heisenberg models are all equiv- 
alent (so fc ^ 25). While we have found within DMRG 
for the 5 = 3/2 TB chain that ttvx is in agreement with 
fc = 3, we find for the S = 3/2 Heisenberg model fc = 1 
(as well as c = 1). 

These findings can be further generalized to SU(A'"), 
SP(A''), or SO(A'')-symmetric WZW models; for instance, 
we expect that ttvx = k{N^ - l)/6 for SU{N)k WZW 
models [107, 108] but leave this topic for future work. 



B. Numerical Results 

1. Zero Temperature 

One can check the previous analytical results in 
large-scale numerical simulations using DMRG or QMC. 
DMRG results have already be presented in Ref. [53] for 
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FIG. 12: Quantum Monte Carlo (QMC) results for the fluc- 
tuations J-'(L/2) versus the total system size L for periodic 
chains of lengths L = 16, 32, 64, 128, 256. Different col- 
ored symbols are for different values of the Ising anisotropy 
of the XXZ chain. For the XX point A = 0, QMC data 
(black circles) are compared to exact diagonalization results 
(red crosses). Logarithmic scaling is observed in the critical 
regime A < 1 whereas for A > 1 TiL/l) tends to saturate at 
large L. Inset; schematic picture for the periodic ring where 
subsystem A is taken with x = L/1 sites. 



the spin-1/2 XXZ chain in the critical regime, and Ap- 
pendix C gives a detailed account of how cumulants of 
charge fluctuations can be determined in DMRG. Here 
we show QMC data obtained at zero and finite tem- 
perature using the Stochastic Series Expansion (SSE) 
method [109, 110]. As in DMRG, it is rather straightfor- 
ward to compute the fluctuations of the magnetization 
of a subsystem A embedded in a spin chain with PBCs. 
We choose the simplest partition, cutting a length-L ring 
into two halves as depicted in Fig. 12, and compute 



T{L/2) 



L/2 



L/2 



(3.35) 



QMC results for the ground-state expectation value 
T{L/2) obtained using T — > extrapolations (see be- 
low) are shown in Fig. 12 for various anisotropics A and 
sizes L = 16, 32, 64, 128, 256. For A = we suc- 
cessfully compared QMC results with exact diagonaliza- 
tion performed using the free-fermion representation. In 
the entire antiferromagnetic critical regime A 6 [0, 1] the 
scaling of J- is logarithmic, while J- tends to saturate 
once A > 1, as clearly visible in Fig. 12. Assuming a 
logarithmic scaling of the form 



^2j-(L/2) =irefflni + /i 



(3.36) 



where /i is a non-universal constant and Kf.s{L) is the 
effective (size-dependent) LL parameter used to fit the 
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FIG. 13: Prefactor Kes of the logarithmic scaling Eq. (-3.36) 
for the fluctuations of the ID XXZ chain computed with 
QMC, see Fig. 12. Kes was extracted by fitting QMC data 
for L e [2P-\ 2P] with L = 128 (p = 7) and L = 256 (p = 8). 
The blue line shows the exact Bethe Ansatz result for A'(A) 
Eq. (3.22). 



QMC data for L € [2p-\ 2^]. Results for K^s as a func- 
tion of A are shown in Fig. 13 for L = 128 (p = 7) 
and L = 256 (p = 8), and are compared with the Bethe 
Ansatz result in Eq. (3.22). There we see three impor- 
tant features: (1) In the gapped regime (A > 1) the log- 
arithmic scaling is of course not present and we find that 
the effective LL parameter rapidly vanishes with L. (2) 
For A < 0.8, the finite-size effects on K^g appear to be 
very small since numerical estimates and Bethe Ansatz 
results almost coincide. This turns out to be natural 
because we know from the renormalization group (RG) 
equations [111] that Kes{L) goes to its infinite value K 
as L"" with exponent [112] a = 8{K - 1/2) which be- 
comes less than 1 if K < 5/8, i.e., A > 0.9. (3) At the 
isotropic Heisenberg point A = 1, a marginal operator 
in the effective field theory induces finite size logarithmic 
corrections to the prefactor Kf.g which makes the con- 
vergence to the thermodynamic limit K = 1/2 slower. 
The same fixed point without such finite-size effects, the 
Haldane-Shastry model, has already been discussed in 
the previous section. 



2. Finite- Temperature QMC Results 

Fluctuations at finite temperature can be compared to 
the Curie constant of the uniform susceptibility. 



where the uniform susceptibility is defined as usual by 

' dm,; 



X{T) 



(3.38) 



h=0 



^-lj2(S!S^)-i2(^t)A ■ (3.39) 
\i,j=i i=i J 



The fluctuations within a subsystem A containing I spins 
reads 



which can be rewritten as 



T{e,T) = £C{T)-J2J2(^!s; 



(3.40) 



(3.41) 



ieAjeB 



The above expression contains useful information. At 
high temperatures, the second term on the righthand 
side of Eq. (3.41) vanishes because there are no corre- 
lations between distant spins. In the limit T — > co, 
therefore, the fluctuations T coincide with £C, as can 
be seen in Fig. 14 where flnite-T QMC results are pre- 
sented for T{L/2,T) and x{T). More precisely, the main 
panel of Fig. 14 shows QMC results for F{L/2,T) and 
L/2 X C(T) obtained for isotropic Heisenberg rings of 
various lengths L = 32, 64, 128, 256. The extensive Curie 
constant L/2 x C{T) starts to deviate from T{L/2,T) 
for T r-^ J. Below this energy scale C goes to zero, as ex- 
pected for an antiferromagnet, and decreases smoothly 
down to very low temperatures where it saturates to its 
ground-state expectation value (3.36) when the tempera- 
tm-e T becomes less than the flnite-size gap A{L) ~ J/L. 
To be more quantitative, let us deflne the "quantum con- 
tribution" 

JFq (£, T) = J-(£, r) - £C{T) (3.42) 
ieAjeB 

which measures the quantum part in the fluctuations 
since thermal contributions are subtracted. Of course 
Tq and coincide at T = 0, as visible in Fig. 14. At 
high temperature, correlations between A and B are very 
short-ranged, being governed essentially by ^(T) ~ l/T. 
In this limit, we therefore expect 



T- 



(3.44) 



CiT) = Tx{T), 



(3.37) 



which is indeed observed numerically as shown in the 
inset of Fig. 14. When T decreases, Tq{L/2, T) increases 
logarithmically: 

J-Q (L/2, T) = -^\n-+ const. (3.45) 

up to saturation when the temperature becomes smaller 
than the flnite size gap A(L) ~ J/L. This will be further 
discussed in the experimental section below. 
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FIG. 14: Finite temperature QMC results for SU(2) Heisen- 
berg rings. Main panel: T{L/2,T) vs. T/J for various 
lengths L = 32, 64, 128, 256 (symbols) plotted together 
with the subsystem Curie constant L/2 x C{T) (dashed lines) 
over a wide temperature range. Inset: Difference Tq{L/2,T) 
from Eq. (3.43) plotted versus T/J. The dashed pink line is 
7 - (lnT)/(27r^) with 7 = 0.15, and the full black line is the 
high-T limit ~ 0.125/T. 



C. Gapped Systems 

Until now we have only considered gapless models 
where we could find universal results for the fluctuations. 
Gapped models are, in contrast, non-universal since a 
new length scale is introduced by the gap. Neverthe- 
less, gapped models possess interesting aspects from an 
entanglement point of view. In gapped one-dimensional 
systems the entanglement entropy obeys a strict area law 
[113], and we expect the same for the fluctuations due to 
the finite correlation length and exponentially decaying 
correlation functions. In one dimension, the scaling be- 
havior of a strict area law is a constant since the region 
at the cut consists of a single point (or two for PBCs). 
For a large enough subsystem A ( "large" being governed 
by the size of the gap) it is therefore sufficient to consider 
the cut at the boundary between both subsystems. 

As a paradigm for gapped systems we consider the 
spin-S" Affleck-Kennedy-Lieb-Tasaki (AKLT) chain [114] 
because we can evaluate all the quantities rigorously. The 
AKLT states, also called valence bond solid (VBS), serve 
as model wave functions for Haldane-gap spin chains 
[115]. Most recently, the AKLT models have attracted 
some interest since they can be seen as topologically non- 
trivial states of matter [116, 117]. The spin S AKLT wave 
function can be most conveniently written as 



m = n 



using Schwinger bosons [118, 119]: |t) = 



(3.46) 
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FIG. 15: (Color online) Fluctuations (diamonds) and entan- 
glement entropy (dots ) for the spin 1 AKLT chain with 
L = 50. Blue (green) data points correspond to PBCs 
(OBCs). The curve fitting the EE for PBCs (blue dots) corre- 
sponds to Eq. (3.48), the black curve fitting the EE for OBCs 
(green dots) is the analytic expression Eq. (3.47). The pur- 
ple line fitting the fiuctuations for PBC (blue diamonds) cor- 
responds to Eq. (3.51). Green data points are obtained for 
OBCs and 5tot = 1, thus providing an example for an AKLT 
state with nonzero magnetization. The bulk values of the 
curves are (from top to bottom) 21n2, ln2, 1/2, and 1/4. 



II) = ftt |o) while the spin 1 states are defined as 
1,1) = l/v/2(at)2 |0), |1,0) = at6t|o), and |1,-1) = 
1/72(6^)2 |o). A singlet bond reads {a\b] - b\a]) |0). If 
we perform the partitioning of the system into A and 
B, the cut inbetween is nothing than an intersected sin- 
glet bond. The entanglement entropy of an intersected 
singlet bond is In 2, the fluctuations are 1/4. Thus we 
might expect the entanglement entropy (fluctuations) to 
saturate to the value In 2 (1/4) for OBC and to the value 
2 In 2 (1/2) for PBC. In the following, we wifl see that for 
OBC this guess is too naive. For clarity we restrict the 
discussion for the moment to the spin 1 case. For PBC 
the AKLT state is a spin singlet, for OBC, however, free 
edge spins (called dangling spms) appear [114]. The dan- 
gling spins with 5=1/2 couple either into a singlet or a 
triplet, yielding a four-fold degenerate ground state. The 
intersected valence bond at the cut contributes always 
In 2 (1/4) to the entanglement entropy (fluctuations). In 
addition, there might be a contribution from the phys- 
ical boundary (i.e., the other edge of subsystem A) if 

state, i.e., either 
For the channels 
either |tt) or 



the dangling spins form a S^^ 

(in> - i;t))/v^ or (in) + i;t))/V2. 



±1, i.e 



with finite magnetization S^^ 
\\,\), the physical edges do not contribute to either en- 
tanglement entropy or fluctuations. 

Exact expressions for the von Neumann entropy have 
already been derived [120-123]. For OBCs the entangle- 
ment entropy is given by 



5"'^^(£) =ln2 



(3.47) 



where ^ ^ = In 3 is the correlation length of the AKLT 
state. This expression, the black curve in Fig. 15, per- 
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The simplest way is to consider the reduced density ma- 
trix, and one will find for the von Neumann entropy 
S = ln(l + S) per cut [122]; similarly, we find that the 
fluctuations are given by T = S{S + 2) /12 per cut in 
agreement with the asymptotic value of (3.51) as well as 
our DMRG findings, see Fig. 16. 

These exact results might be extended for various gen- 
erahzed AKLT states [119, 125, 126] where the SU(2) 
symmetry of the electron spin is replaced by SU(n), 
SP(n), or SO(n) symmetry. 
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FIG. 16: (Color online) Fluctuations (blue dots, bottom) and 
entanglement entropy (red points, top) for the spin 2 AKLT 
chain. The purple line corresponds to Eq. (.3.51), all data 
points obtained within DMRG {L = 50 and PBC). J- satu- 
rates to 4/3 and 5 to 2 In 3. 



IV. CONCLUDING REMARKS ON BIPARTITE 
FLUCTUATIONS 

A. Comparison to Several Proposed Measures of 
Many-Body Entanglement 



fectly matches the DMRG data (green dots). For PBCs, 
the von Neumann entropy is given by [121, 122] 

S^^^(£) = -SXa'^-rXa- Xs^nXB , (3.48) 

where = i^^S^ = ^^'t'^tr"' ; 

L is the system size and p = —1/3. Again, the agreement 
between this anaytical expression and the DMRG data 
is excellent (see top curve in Fig. 15). In what follows 
we compute the fluctuations for general spin S. We start 
with the static spin-spin correlation function [118, 124] 

{StSp - (Sf)(^|) = i^±il!(-l)-^e-l-.'-|/« (3.49) 

with correlation length = In (1 -|- 2/5). The sum 
over i,j can be carried out by writing 

-^aklt(^) = ^^^±^£ + 2j2ie-k)(-e-'^^)' , (3.50) 

k=i 

which is a geometric series and therefore easily summed: 

-Faklt(0 = ^^^^[1 - (-l)'e-^/^]. (3.51) 

The last term in brackets vanishes for ^ — > cx:i, so the 
fluctuations reach a constant. In Fig. 15, J^aklt(^) is 
plotted together with the DMRG data. We still find a 
constant ratio S/T = 4 In 2 for the spin-1 VBS state. 
The ratio depends purely on the cuts (intersected valence 
bonds) rather than on microscopic details. Remarkably, 
the ratio is the same for PBCs and OBCs as it is the case 
for conformally invariant models. 

For the generalization to higher-spin AKLT chains one 
deals with S intersected singlet bonds per cut rather than 
a single one. The bonds are, however, not independent 
and it is not simply a factor S which renders the cor- 
rect entanglement entropy or fluctuations, respectively. 



As noted in the Introduction, the difficulty of com- 
puting the entanglement entropy (1.1) for most many- 
body systems described by microscopic Hamiltonians has 
led to the search for alternative measures of entangle- 
ment which behave similarly to the entanglement en- 
tropy. Most of the measures behave quite similarly in 
one dimension to the entanglement entropy and to each 
other, while in higher dimensions they generally differ. 

One such quantity is the valence bond (VB) entangle- 
ment entropy for SU(2) quantum spin systems [127, 128], 
which generalizes the idea, first used in the context of 
random spin chains [74, 86, 87], that for a pure valence 
bond state whose wave function is a product of disjoint 
singlet bonds between two spin-l/2s the entanglement 
entropy of a subsystem is essentially the number of sin- 
glet bonds that cross the boundary between the subsys- 
tem and the remainder of the system, with each singlet 
l^'s) = (ITa^b) - \iAtB))/V2 contributing ln2 to the 
entanglement entropy. For more general wave functions 
which are superpositions of pure VB states the number 
of such crossings can be averaged to obtain a well-defined 
(with certain restrictions [127]) measure of entanglement. 
In one-dimensional critical systems the scaling of the VB 
entanglement entropy has been shown to be log- 
arithmic with subsystem size, though without the cen- 
tral charge-dependence of the von Neumann entangle- 
ment entropy [129]. It is also now clear that in higher 
dimensions the VB entanglement entropy does not nec- 
essarily scale exactly as the von Neumann entanglement 
entropy; for the spin- 1/2 Heisenberg antiferromagnet on 
the square lattice in two dimensions, for example, the 
von Neumann entanglement entropy obeys a strict area 
law with an additive logarithmic correction while the VB 
entanglement entropy contains a multiplicative logarith- 
mic correction [45, 130]. Interestingly, however, the strict 
area law can be recovered if a modified definition of en- 
tropy based on loops of the transition graph is used in 
the valence-bond basis [131]. 

For pure VB states the fluctuations can also be un- 
derstood in a similar manner: for each singlet we can 
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consider the spin fluctuation Ta = {{S\f) = 1/4.. Then 
for a pure VB state with n bonds crossing the bound- 
ary, S = = (4 In 2)7" = nln2. This idea is explored 
further below. 

Another suggested measure of entanglement is the log- 
arithmic bipartite fidelity [132] F = -\n{\(AU B\A (g) 
B)!^), where the state \AU B) is the ground state of a 
many-body Hamiltonian and \A(g B) = \ A) \B) is the 
ground state of the same Hamiltonian with all interac- 
tions between the two subsystems removed. The argu- 
ment of the logarithm, |(^ U B\A(Ej B)\^, has an inter- 
pretation as the probability of finding the ground state 
energy if the energy is measured just after connecting 
two previously disconnected subsystems. What makes 
the logarithmic bipartite fidelity particularly interesting 
is the fact that not only is the scaling with system size 
logarithmic in one-dimensional gapless systems but the 
prefactor is proportional to the central charge of the un- 
derlying CFT just as for the entanglement entropy. Thus 
the appearance of the central charge is not unique to the 
entanglement entropy. On the other hand, unlike the 
entanglement entropy the finite-size result cannot be ob- 
tained from the infinite result by a simple substitution 
arising from a conformal transformation, and the quan- 
tity remains difficult to both measure experimentally and 
compute in numerical simulations beyond one dimension. 



stead: 



aa; In a; + bx 
b'x 



(Neel), 
(VBS). 



(4.2) 



If X denotes the perimeter of the boundary between 
the two subsystems, QMC estimates are a ~ OA and 
b ~ 0.02 [127] for the Heisenberg antiferromagnet on the 
square lattice. 

For a single VB state \ip), each virtual bond b = [ij] 
connecting two lattice sites i and j can be assigned by 
its occupation number rif, such that 

^ ^ J 1 if 6 is occupied by a singlet, 
1 otherwise. 

Then the VB entanglement entropy for such a state \ip) 
between two subsystems A and B is defined as follows 



5^B(^) = ln2x^nb(^), 



(4.4) 



where the sum is performed over all the bonds 6a con- 
necting the two subsystems A and B. Of course, only 
the occupied bonds will contribute 1. This obviously 
puts an upper bound on iS^^ for a subsystem of A'^ sites: 



5„X£, = iVln2. 



2. Links with T 



B. Relation to the Valence Bond Entanglement 
Entropy 

1. Basic Notions 

The VB entanglement entropy was introduced as an 
alternative measure of bipartite entanglement in SU(2) 
symmetric quantum spin systems [127, 128, 133]. While 
conceptually distinct from von Neumann and Renyi en- 
tropies [43, 45], VB entanglement entropy (and alterna- 
tive definitions in the VB basis [131]) has been shown to 
be an efficient entanglement witness and a powerful tool 
to detect quantum phase transitions [127, 134] in SU(2) 
spin systems. Among several interesting features, the VB 
entanglement entropy was computed analytically for the 
critical spin-1/2 Heisenberg chain in Ref. [129] by Jacob- 
sen and Saleur, who found a logarithmic scaling of the 
form 



41n2 . 



ln.T + ri, 



(4.1) 



where ri is a constant. Such a scaling turns out to 
be accidentally close to the CFT result for the entan- 
glement entropy [7] where the prefactor 1/3 is similar 
to the 4(ln2)/7r^ ~ 0.28 above. In two dimensions, 
the Neel-ordered spin-1/2 antiferromagnetic Heisenberg 
model was shown to exhibit a log-modifed area law [45, 
127, 128]. If an enforced bond dimerization leads to a 
gapped VB solid state, a strict area law was found in- 



For a SU(2)-symmetric system, if the singlet ground 
state is a unique dimer covering |$o) = Iv'); the two-point 
correlation is simply given by 



($o|5f^||$o)=-^n[,,]. 



(4.5) 



Using the fact that at zero temperature (see Eq. 3.41) 
ieAjeB 



we have the simple result 



T = 



S 



VB 



41n2 



(4.7) 



In the random singlet state where S^^{x) = S"^^(x) = 
[(ln2)/3] lux, it yields T{x) = (1/12) Inx. 

However, a general singlet state can be a complicated 
superposition involving a large number of dimer cover- 
ings, such that one can expand any singlet state on a 
(overcomplete) valence bond basis 



l*o) =X!Ap|<^p). 



(4.8) 



In such a case the VB entanglement entropy has been 
defined in Refs. [127, 128] as 



5^^($o) = 



(4.9) 
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From such a definition, one sees immediately that the 
equality (4.7) is not valid anymore. Indeed, if one defines 



ieAjeB 

from Eq. (4.9) we get 



J2p K 



whereas 



p p' 



(4.10) 



(4.11) 



(4.12) 



Clearly we see that, except for trivial general su- 

perposition of the form Eq. (4.8) yields F ^ S'^^/iA In 2). 
An illustration based on a simple 4-site example is pro- 
vided in Appendix D. 

Despite their fundamental differences, it is quite inter- 
esting to notice that for the Heisenberg chain, when com- 
paring the exact expression of Jacobsen and Saleur [129] 
Eq. (4.1) for S^^ with F, one sees that the leading log- 
arithmic term of 



-Fidhaf(^) = ^ln£ + /i. 



(4.13) 



is directly proportional with a factor 8 In 2 (to be com- 
pared with 4 In 2 for a pure valence bond state). However, 
this coefficient does not hold when moving away from the 
Heisenberg point along the XXZ critical line [129]. 



C. Non-Equivalence of Fluctuations and 
Entanglement Entropy 

As noted in the Introduction, the equivalence of the 
charge statistics with the entanglement entropy for non- 
interacting fermions, while suggestive, is particular to the 
nature of non-interacting particles: here entanglement 
entropy arises solely from the uncertainty of the posi- 
tion of particles and therefore particle-number fluctua- 
tions within a given volume encodes the full information 
needed to determine the entanglement entropy. We have 
shown that many essential features are retained in the 
fluctuations for interacting systems, especially when the 
theory is essentially gaussian as in one-dimensional Lut- 
tinger Liquids. It is worth noting, however, that, even 
in principle, it is unlikely that there exists a universal 
series similar to Eq. (2.1) for interacting systems with 
more complicated coeflicients. This can be seen by ob- 
serving the simple fact that when calculating the fluc- 
tuations (or cumulants) from the reduced density ma- 
trix, one essentially averages over the eigenvalues in each 
particle-number sector, so that the individual eigenval- 
ues cannot be discerned. Stated another way, the full 
counting statistics cannot, in general, contain the same 



amount of information as the entanglement entropy. The 
example for DMRG in Appendix C makes this clear. 

To see this more explicitly, it is already known that 
for an ideal Bose gas the relation between entanglement 
entropy and fluctuations is actually S ~ In [39, 135]. 
In the absence of particle number conservation in the 
full system, moreover, the arguments given in the Intro- 
duction for the symmetry between subsystems no longer 
holds and the fluctuations are not expected to behave as 
the entanglement entropy. For the quantum Ising chain 
in transverse fleld, for example, the fluctuations of do 
not scale logarithmically [73]. Of course, for this partic- 
ular case one must instead relate the entanglement en- 
tropy to fluctuations of the appropriate quasi-partides 
(see Supplementary Material). 



D. Accessible Entropy and Charge Fluctuations 

Charge conservation restricts the possible local oper- 
ations available when using fermionic modes for quan- 
tum information processes purposes. A modification of 
the usual entanglement entropy, which quantifies the use- 
ful, accessible entanglement was proposed in Ref. [136]. 
There, it was suggested to quantify the accessible entan- 
glement entropy S'^^ by averaging entanglement entropy 
over super-selection sectors (see also Refs. [79, 137, 138]). 
In Ref. [139] it was suggested that constraints on local 
operations may be treated as a resource for hiding in- 
formation in correlations which are blocked from local 
probing. 

In Ref. [140] , charge fluctuations have been used as an 
estimate of the difference between accessible and inacces- 
sible entropy via the inequality: 



Sa-AS< ST < Sa, 



AS = i log 



27re Ca + 



where 



p„.™ = Tr(n;^$§n^pn;^®n^), 



(4.14) 
(4.15) 

(4.16) 
(4.17) 



where H^ project onto sectors with flxed particle number 
n in A (i.e., on states ip in such that = nij}). 

Similarly, H^ projects on sectors with A''^ = m m B. 

An immediate consequence of the examples we have 
considered in the previous sections is that the C2 scal- 
ing can only lead to sub-leading logarithmic (or log log 
in one dimension) corrections to the accessible entangle- 
ment entropy. 



E. Experiments 

In principle, bipartite fluctuations are accessible in ex- 
periment, which was part of the motivation for studying 
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them in place of entanglement entropy, which remains 
the quintessential measure of "quantum-ness" in a many- 
body system. The latter, however, has always suffered 
from lack of measurability (but see Ref. [42]) despite the 
great theoretical interest stemming from fundamental re- 
sults. In this work we have shown how for non-interacting 
fermions the fluctuations are equivalent to entanglement 
entropy [36, 40], while for interacting systems we trade 
non-equivalence for measurability. It remains to show 
how this can be accomplished in real systems, and below 
we propose and discuss three different systems where T 
might be measured. 

1. Quantum Point Contacts 



superconducting piece 




A 



- — Quantum magnet 

top view side view 



As discussed in Sec. II CI, a quantum point contact 
(QPC) is a beam splitter with tunable transmission and 
reflection that serves as a "door" between electron reser- 
voirs [36, 40]. The current fluctuations J^{t) measured 
at this QPC are exactly our bipartite fluctuations when 
replacing the temporal window t by the spatial extent x. 
Due to space-time duality, such a replacment is always 
allowed for conformally invariant sytems [36] , even when 
they are interacting. Although the two situations are not 
identical, the QPC measurements being more closely re- 
lated to the concept of local quantum quenches [75, 141], 
QPCs with free-electron reservoirs offer the opportunity 
to unambiguously measure entanglement entropy using 
the formalism of Sec. II. The flrst few cumulants, in par- 
ticular C2, have already been successfully measured [58]. 

2. Cold Atoms 

The second system where the fluctuations are avail- 
able are cold atomic quantum gases loaded into optical 
lattices. To be more speciflc, we consider the recent ex- 
periments where site-resolved images of atoms in opical 
lattices have been taken for the first time with a quantum 
gas microscope [57, 61, 62, 142]. In these experiments, 
it was possible to measure the mean and variance (i.e., 
the fiuctuations) of the particle number per site. Strictly 
speaking, current technology allows for the measurement 
of parity number, which corresponds to particle num- 
ber only for free fermions and hard-core bosons (in the 
Tonks regime). It is possible, however, that in future ex- 
periments actual particle number at the subsystem level 
may be measured. 

3. Quantum Antiferromagnets 

Third, bipartite fluctuations can also be measured in 
quantum magnetic systems, with the ^-component of 
spin (or the quantization axis, more generally) play- 
ing the role of particle number in QPCs and cold atom 
systems [143]. We note that Eq. (3.23) essentially defines 



FIG. 17: Experimental setup proposed to extract the fluctua- 
tion Ta of tlie magnetization within a subregion A. A super- 
conducting device (black disk) is placed on top of the quantum 
antiferromagnet (on region B) so that the Meissner effect will 
induce hs = in the antiferromagnet whereas Ha = H. The 
superconducting device is also placed on the other side such 
that only the A contribution of the field-induced magnetiza- 
tion will be measured. 

the Curie constant of the partial uniform susceptibility, 
J-A = T X d{S\) /dhj^lfij^^o where hji is a small external 
uniform magnetic field applied only to region A. Such 
a setup could be realized by applying the magnetic field 
over the entire sample and covering the region B (the 
remainder of the system) with superconducting Meissner 
screens [130]. Such screens would eliminate the external 
field as well as the magnetic response outside the region 
A. By varying the size of A, the scaling of Ta could be 
measured. A setup is sketched in Fig. 17 

F. Conclusion 

We have presented a detailed study of the behavior of 
bipartite fluctuations of conserved U(l) charges such as 
particle number A'^ and spin in various quantum many- 
body systems and proposed that the study of such fiuc- 
tuations is a useful tool for characterizing such systems. 
For systems that can be mapped to a problem of non- 
interacting fermions we have shown that the fiuctuations 
and higher-order cumulants encode all of the information 
needed to determine the commonly used measures of en- 
tanglement, the von Neumann entanglement entropy and 
Renyi entanglement entropies. As we have shown, the full 
entanglement spectrum can also be determined in this 
manner, demonstrating for the first time that, in princi- 
ple, the entanglement spectrum is an accessible quantity. 
In other systems the behaviors of the fluctuations and 
entanglement entropy differ, but the fluctuations still re- 
veal important features of the system while having the 
advantage of being tractable analytically in many cases, 
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computable numerically in both DMRG and QMC, and 
in principle experimentally accessible. We therefore hope 
that this work will motivate further work into under- 
standing more fully the behavior of the fluctuations as 
well as their relation to the entanglement entropy in in- 
teracting systems in one and higher dimensions. 
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APPENDIX A: DERIVATION OF THE RENYI 
ENTANGLEMENT ENTROPY SERIES 

In this appendix we derive Eq. (2.3) following similar 
steps to the derivation for the von Neumann entangle- 
ment entropy in Sec. II A. As shown in the Supplemen- 
tary Material, the Renyi entropies can be written in the 
form 

Sn = ^Tr{ln[M" + (1 - A/)"]}. (Al) 
We first expand the logarithm in Eq. (Al) as 
ln[A/'^ + (1 - Mr] ^ - ± ii^-Mn- (l-Mr]\ 

.7 = 1 

(A2) 

By repeated use of the binomial theorem, we obtain 



Tr{ln[Af" + (l-Af)"]} = -^^ ^^(-l)'"+'-+^i( Ml-' I ""')Tr(A^"'^+=). (A3) 

j = l m=0 r=0 s=0 V^J \ ^ / \ / 

Before using Eq. (2.16) to write this in terms of the factorial cumulants, we must show that Fq does not contribute, 
i.e., that the expression vanishes when nr + s = 0. Since n > 1 and the sum over m vanishes when r = s = this 
is indeed the case. Then we can substitute Eq. (2.16), and moreover use relation between factorial and ordinary 
cumulants (see Supplementary Material) to get 

oo j j-m nm nr+s \/ \ o / i ^ 

( J \ I J ~ ''^\ j nm^^ bi(nr + s,w) ^ 

j = l m=0 r=0 s=0 w=0 



Tr{ln[Af" + (1 - A/)^} ^ - ^ E E E E ("D'""^"^" 7 L ~ 7 Zll-iy '^- ^^^^ 



It is now possible to introduce a cutoff R in the outermost sum and write Tr{ln[A/" + (1 - A/)"]} = Pk{n, R)Ck, 

where, after some algebra, 

1 — n r \ r J \m J \nr — sj (s — 1)! 



The coefiicients l3k{n,R) actually vanish for k > nR, so 
we finally arrive at Eq. (2.4). As a check on the alge- 
bra, it is again found that for _R — > oo the coefficients 
Pk {n, R) approach the limiting values derived in the Sup- 
plementary Material, although this time only a numerical 
demonstration is possible. 



APPENDIX B: CONVERGENCE OF THE 
ENTANGLEMENT ENTROPY SERIES 

In this appendix we prove the convergence of se- 
ries (2.1) for the entanglement entropy, which is relatively 
simple to show under the assumption that Al is finite. In 
the derivation of Eq. (2.23) the cutoff K originates from 
the expansion of Eq. (2.9) and hence convergence in the 
limit — > oo corresponds to convergence of the partial 
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sums of Eq. (2.12), i.e., that Eq. (2.12) converges to the 
correct sum for the M of interest. Now, let us first con- 
sider the case where M is a single number p such that 
< p < 1. This is the property satisfied by the eigen- 
values of M for general i\/ for fermionic systems due to 
Fermi-Dirac statistics (strictly speaking, < p < 1, but 
at the extremes we understand the entanglement entropy 
to be zero). The series 



(Bl) 



is absolutely convegent by the ratio test, since 
luiifc^oo |afc+i/afc| is 



lim 



[p(l - pf+^ + p''+\l - p)]/{k + 1) 
\p{l-pf +p''{l-p)]/k 

l-p if0<p<l/2, 
p ifl/2<p<l. 



(B2) 



In both cases the limit is less than 1 and the series there- 
fore converges. Moreover, each term is positive and by 
necessity decreasing so that the partial sums approach 
its limit from below. 

If Eq. (2.9) is evaluated in tlie eigenbasis of M then 
we can apply the above argument to each eigenvalue, 
but in general M is a matrix with off-diagonal elements. 
From the theory of matrices, however, it is known that 
Eq. (Bl) converges for p replaced by a general matrix M 
if < a{M) < 1 where a{M) is the spectral radius or 
maximum eigenvalue of M. Therefore the matrix series 
converges as well. 

We can similarly show that Eq. (2.3) for the Renyi 
entropies converges. In this case the cutoff R arises from 
the series expansion (A2). Again let us assume first that 
M consists of the single entry p, with < p < 1. We 
have 



lim 



|(l-p")-(l-p)n <1 



(B3) 



where p" + (1 — p)" is positive and p" + (1 — p)" < 
p + (1 — p) = 1 for n > 1, so that the series is again ab- 
solutely convergent. Application of the matrix spectral 
radius to M then shows that Eq. (2.3) also converges, 
and recalling the factor of 1 — n left out in Eq. (A2) we 
see that each term of the series for the Renyi entropies 
is positive and decreasing so that increasing K improves 
the approximation from below as for the entanglement 
entropy. 



APPENDIX C: CUMULANTS IN DMRG 

In this appendix we describe how to compute the cu- 
mulants of a conserved U(l) charge within DMRG; we 
show that the fluctuations and higher-order cumulants 
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FIG. 18: Eigenvalues A; of the reduced density matrix for a 
subsystem A, obtained within DMRG. The data is produced 
for an XXZ spin chain, Eq. (3.20), with A = 0.2 and L = 32. 
(Left) Negative logarithm of the eigenvalues, — In A;, i.e., the 
real-space entanglement spectrum. (Right) The same data vs. 
the spin quantum number for subsystem A. As expected, 
the SU(2) multiplets are not present due to broken SU(2) 
symmetry for A = 0.2. 



can be computed with no additional computational ef- 
fort relative to existing implementations. One can of 
course do this by computing all possible two-point cor- 
relators {xiXj) — {xi){xj), where x = or x = n, re- 
spectively, inside subsystem A; by performing the sum 
= "l^i jeAi^i^j) — {xi){xj) one obtains the fluctua- 
tions in A. It is, however, costly to compute the nec- 
essary correlators within DMRG when the system size 
becomes large. We therefore present a much more el- 
egant way to get the fluctuations directly without any 
need for computing the correlators: the fluctuations, and 
indeed all higher-order cumulants, can be extracted sim- 
ply from the reduced density matrix of the subsystem 
in the usual DMRG implementation with no additional 
computational effort. 

In each DMRG step the full system is split into two 
blocks, A and B, with the reduced density matrices pA 
and pb, respectively. Without loss of generality, we can 
use the eigenvalues of pa to compute the von Neumann 
entropy in the Schmidt basis. 



Sa = Ailn Aj, 



(01) 



so that the entanglement entropy is a natural byproduct 
of DMRG, since all A^'s are computed in every DMRG 
step. 

As mentioned in the Introduction, the reduced density 
matrix is usually stored in block-diagonal format accord- 
ing to the U(l) quantum number (i.e., either the particle 
number n or the spin projection Sz) [38]. In other words, 
the eigenvalues A^ are in fact labeled and we should write 
a'^' where x is the U(l) quantum number. Now we can 
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write the fluctuations in the Schmidt-basis as 



(C2) 



where the sum is restricted to a fixed value of the 
quantum number x. Higher moments are computed anal- 
ogously while the cumulants are determined from the mo- 
ments via well-known formulas. 

The reader should notice that we refer here to x (e.g., 
Sz) of block A (or B) but not to the total x of the whole 
system. To avoid confusion one could write x^ instead of 
X. Eq. (C2) implemented in the DMRG code provides an 
efficient way to compute the flucutations. We present a 
representative example in Fig. 18 for the eigenvalues of 
an XXZ spin chain, Eq. (3.20) with A = 0.2 and L = 32 
where the bipartition is made at L/2. It is clear that 
the fluctuations arise quite naturally in DMRG compu- 
tations. 



APPENDIX D: ILLUSTRATION OF THE 
DIFFERENCE BETWEEN T AND 5^^ 

In order to illustrate the fundamental diff'erences for a 
SU(2) symmetric model between the fluctuation J- and 



the VB entanglement entropy , let us take this simple 
normalized 4-site singlet state 



1 
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(Dl) 



Here a dimer — connects two sites belonging to differ- 
ent sublattices A and B, and is written with the following 
convention | — = | tAj-B - 4-AtB)/v^- 

The bipartite entanglement shared across a vertical 
boundary drawn in the middle of the 4 site system is 
then easy to compute. Regarding VB entanglement en- 
tropy, I ) contributes x In 2 and | ) 
contributes 2 x In 2. Using the deflnition Eq. (4.9), we 
simply get 



0-h21n(2)/x/3 



= In 2, 



whereas computing the fluctuation yields 

cVB 

T=- 

61n2' 

which is different from Eq. (4.7). 



(D2) 
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I. THE REDUCED DENSITY MATRIX OF NON-INTERACTING FERMIONS IN 

EQUILIBRIUM 



In this section we use the special form of the reduced density matrix for non-interacting fermions in 
equilibrium to derive the following essential formulas for the von Neumann entanglement entropy S, 
Renyi entanglement entropies Sa, and the generating function x(-^) for particle number fluctuations 



in a system of non-interacting fermions: 



S = -Tr[Mn InMn + (1 - Mn) ln(l - Mn)], (1.1) 

= -^Trln[A/^ + (1 - Mn)% (1-2) 
I — a 

x(A) = det[l + (e^^-l)A/n]. (1.3) 

As explained below, Mq is a Hermitian matrix with dimension equal to the number of fermionic modes 
in the subsystem O of interest. Eqs. (1.1)-(1.3) apply to any system of non-interacting fermions, but 
are particularly useful when applied to a subsystem of a system of non-interacting fermions. All of 
the results in this section are known in some form, especially from Ref . [1] , but we provide a detailed 
presentation for clarity and to point out relevant features of the derivation along the way. 



A. Without Pairing 

We begin with the simpler case where fermion number is conserved. Non-interacting fermions in 
any dimension are described by a Hamiltonian of the form 

^ = ^alf/,,a„ (1.4) 

where Ui is the fermionic annihilation operator for degree of freedom i (usually i labels a site on the 
lattice, but more generally may carry a spin index) and Hij is a Hermitian matrix. The Hamiltonian 
can be diagonalized by a unitary matrix U whose columns are the eigenvectors of the matrix H such 
that A = HU is diagonal with matrix elements A^^k' = ^k^kk'- Then the Hamiltonian can be 
written in canonical form as 

H = Y,^khlhk, (1.5) 

k 

bk = J2ulA- (1-6) 
At equilibrium the density matrix of the system is 



(1.7) 



Tr e-^H 

and in particular the correlation matrix, or Green's function, at temperature T = 1//3 is given by 

= (a]a,) = Tr(pa]o,) = ^ U.kfmkWl, (1.8) 

k 

where 

m - ^ (1-9) 

is the Fermi occupation function. More simply, we may write 

M = fm). (1.10) 

Clearly A/ is hermitian. Although we will almost always be concerned with the zero-temperature 
ground state, i.e., f{/3Ak) = 6{—Ak), in this section we will keep the discussion general. This will 
also allow us to briefly explore the behavior of the entanglement entropy and fluctuations at finite 
temperature in the main text. 
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Suppose we are interested in a subsystem fl of the total system, usually representing a spatial 
subregion of the lattice. All properties of this subsystem are determined by the reduced density 
matrix pn, which may be found using the fact that the reduced density matrix can be written in the 
thermal form [1] 

Pn = ^ ^ , (1.11) 
Hn=J2 ^UHnh&j, (1.12) 

where the "temperature" Tq = l//?o of the subsystem can be set to 1 without loss of generality. 
We note that Eq. (1.11) by itself, without Eq. (1.12), is completely general and applies to any sys- 
tem, whether composed of non-interacting fermions or not: it is simply a way to write the positive 
semi-definite Hermitian density matrix, and the Hamiltonian Hq is known in this context as the 
entanglement Hamiltonian. The eigenvalues of Hq, called the entanglement spectrum, have recently 
been of considerable interest in characterizing topologically-ordered phases [2]. The additional con- 
straint expressed in Eq. (1.12), i.e., that the entanglement Hamiltonian for non-interacting fermions 
also describes a system of non-interacting fermions, is what makes the following computation possible. 
Briefly, this is a consequence of Wick's Theorem: since the system is composed of non-interacting 
fermions and Wick's Theorem must apply for all correlation functions, in particular to correlation 
functions within subsystem f2, the reduced density matrix for subsystem fl must be of the above 
form [1]. 

Since the density matrix described in Eqs. (1.11), (1.12) is identical in form to the original problem, 
we see that if u plays the role of U, A the role of A, and Cq the role of bk in the previous calculation, 
the correlation matrix in subsystem fl is given by 



- '"^"-^ -Ti{pna]a,) = J2urkf{\k)ul^, i,jefl. (1.13) 



fc 



As for M, we have that m = f{h). Now, it is clear that m^- must agree with Afy for i,j in subsystem 
fl. Moreover, the eigenvalues of m are precisely f{\g), so that if Vg are the eigenvalues of the matrix 
Mq obtained by restricting (or projecting) AI to A, then 

fi\) = (1.14) 

With Eq. (1.14) we can now derive Eqs. (1.1)-(1.3). First, the von Neumann entanglement entropy 
S = — Tr(pn In pn) is simply the ordinary "thermodynamic" entropy of a system of free fermions with 
energies A, at inverse temperature /? = 1, which implies 

'S = ^ff2(/(A,)) = ^ff2(i^,), (1.15) 

1 Q 

where 

H2ix) = -xlnx - [l- x)\n(l- x) (1.16) 

is the binary entropy fimction. Eq. (1.1) is the matrix form of Eq. (1.1.5). 

The Renyi entropies are usually not considered in thermodynamic applications so we carry out 
a more explicit computation by writing out the density matrix (1.11) in the basis in which the 
Hamiltonian is diagonal: 

XX _|_ g A, 

In the g-occupation basis |0q), \ lg) we have explicitly 

1 - /(A,) \_<>,fl-'^g 
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so that 



1^(Pn) = IlTr (^^ o'^'^ Ho) = IlK + (1 - ^9)"]- (1-19) 

Therefore 

5a = ^ ^ ln[< + (1 - i^X] (1-20) 

9 

whose matrix form is Eq. (1.2). We can also confirm that Eq. (1.20) agrees with Eq. (1.15) in the 
limit Q — > 1. 

The statistics of the number operator in subsystem fl, 

7Vn = ^ala, = ^cte„ (1.21) 

may be computed in a similar fashion. We are interested in the generating function 

e(^) = (2*")=Tr(po>), (1.22) 

from which the standard generating function xW is obtained by letting z = e'^. The factorial 
generating function Xf{^) is similarly given by letting z = X + 1. Using the explicit form of the 
reduced density matrix in Eq. (1.18) we have 

g(^) = Trn 1 + (1-23) 

9 

9 ^ ^ 

= l[[l + {z-l),y,], (1.25) 

9 

which in matrix form is Eq. (1.3) with background charge Q = and the substitution z = e*^. 

It is worth emphasizing that, given the assumptions in Eqs. (1.11), (1.12), the only additional input 
to the previous computation is the correlation matix of the system restricted to subsystem Q, i.e., 
Mq. This then determined the fermionic energies Xq in the subsystem, which in turn determine 
the entanglement entropies and number fluctuations. It is not necessary, in particular, that the 
correlation function arise from a system of fermions in equilibrium; the procedure for deriving the 
relevant equations are completely general as long as Wick's theorem applies. 

For example, this allows us to apply the same formulas to non-equilibrium problems where the 
correlation matrix M evolves in time, usually described by a unitary time evolution operator in the 
space of single-particle modes. The case of a quantum point contact (QPC) was originally described 
in Ref. [3], though as we shall see, this description used a flawed formula for relating the entanglement 
entropy to the cumulants of charge transfer. 



B. With Pairing 



For completeness we derive Eqs. (1.11), (1.12) in the presence of pairing between fermions, which of 
course arises in applications involving superconductivity. In this case particle number is not conserved 
so the interpretation of "number fluctuations" is somewhat problematic, especially for experiments. 
It does, however, clarify the source of entanglement entropy in such systems so that we present 
the derivation. Pairing is taken into account by considering the more general quadratic fermionic 
Hamiltonian 



H = ^ alA^jaJ + ^(a|B,ja] + h.( 



(1.26) 
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where A is a hermitian matrix and B is an anti-symmetric matrix. For simplicity we assume that 
both A and B are real. To diagonalize the Hamiltonian we perform the transformation [4] 



(1.27) 



Letting 



(l) = U + V, 
^ = U -V, 



the transformation (1.27) is both canonical and diagonalizes the Hamiltonian as 

H = Y,Kkblb^+H^, ifo = |Tr(A-A), 



(1.28) 
(1.29) 



(1.30) 



with h.kp = ^kSkp, if the following coupled set of equations are satisfied by the orthonormal matrices 
(f) and ■i/': 



{A + B)cj, = iljA, 
{A - B)iP = 0A. 

We define the generalized equilibrium correlation matrix 



G^J = {{ai - aj)(aj + o])) = ^ T/'ifc ( tanh 



(1.31) 
(1.32) 



(1.33) 



I.e., 



Then 



G = V'(^tanh^j<^^ 



tanh' 



/3A 



(1.34) 



(1.35) 



As in the previous calculation, by Wick's Theorem the reduced density matrix for a subsystem fl is 
the same form as the original Hamiltonian: 



-l3nHa 



Pn 



Tr e-z^nfl^n 

Hn = Yl «l(^n)ya,, + li&liBnhja] + h.c.) 



(1.36) 
(1.37) 



where we can again set Tfj = 1. The diagonalization procedure is exactly the same for the Hamiltonian 
of the subsystem as for the original system. Thus if we let (G^G)o be the restriction of G^G to the 
degrees of freedom in f2, then the energy levels A, of the entanglement Hamiltonian Hq corresponding 
to the energy levels Ag of the original Hamiltonian H can be obtained as the eigenvalues of {G'^G)n 



„2=tanh2^. 

g 2 ' 



(1.38) 



(1.39) 
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Comparing to the case without pairing, we obtain all of the same results for the von Neumann and 
Renyi entanglement entropies by setting 



1 - [[G^G)nY'^ 



(1.40) 



This is of course consistent with the fact that G = 1 — 2M in Eq. (1.33) when B = 0, the case without 
pairing. 

One important difference between the case with pairing and the case without, of course, is that 
the fluctuations of particle number in the subsystem do not represent fluctuations of the number of 
original particles. The number operator of quasi-particles in the subsystem Q, is Nq = cjc,, where 
the fermionic operators Cg arise from diagonalizing Eq. (1.37). This is the operator being counted by 
Eq. (1.22). In contrast, the number operator for the original particles is = X^igo o-lo-i, and in the 
presence of a pairing term in the Hamiltonian A^^ ^ A'^^ . 



II. THE COEFFICIENTS OF THE ENTANGLEMENT ENTROPY CUMULANTS 
SERIES IN THE INFINITE-CUTOFF LIMIT 

In this section we perform a calculation that generalizes the computation of Ref. [3] to the Renyi 
entropies and therefore provides, as a limiting case, the von Neumann entropy as well. As noted 
in the main text, the series resulting from such a derivation is not convergent and therefore cannot 
be used to approximate the entanglement entropy with a flnite number of cumulants. The result, 
however, does provide an important check on the values of the coefficients of the series when the cutoff 
is taken to infinity without taking into account the number of cumulants. The derivation, moreover, 
will allow us to point out where this method breaks down. 



1 . Spectral Density Function of Mn 
Let fi{z) be the spectral density function of the correlation matrix Mq, 



fi{z) = Ti[S{z-Mn)] 
1 



= Trj — lim Im — 

\ 7r e^o+ z — AlQ — le 

= — Im lim dzTT\hi{z — Mq — ie)]. 

We recall that has eigenvalues in z 6 [0, 1]. From Eq. (1.3) we have 

Inx(A) = Tr{ln[l + (e^^ - l)Mn]} - (Tr Q)iX, 
so that by factoring out (1 — e*^) in the appropriate manner we can write 

1 



lnx(A) = Tr 
With the substitution 



In 



1 - e' 



Mn 



+ (dim A/n)ln(l 



1 



(Tr Q)i\. 



1 - 

we therefore obtain 

InxiMz)) = Tr[ln(2 - Mn)] - (dim Mn - Tr Q) Inz - (Tr Q) ln{z - 1). 



(2.1) 
(2.2) 

(2.3) 
(2.4) 
(2.5) 

(2.6) 
(2.7) 
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Compaing to Eq. (2.3), we see that the spectral density function can be written as 



ti{z) = -Im[a, lnx(A(2 - iO+))] + (dim Mn - Tr Q)5{z) + (Tr Q)6{z - 1), 



(2.8) 



where we have used 0^ as shorthand for taking the Hmit of e — > In inverting Eq. (2.6) to get A 
in terms of z, we have chosen 



-TT — i hi ( 1 



(2.9) 



Although Eq. (2.6) defines A only within multiples of 27r, this turns out to be the correct phase 
convention. 

As an example, consider the simple case where Mq_ = 1/2 and Q = 0. The spectral density function 
in this case is clearly 



We have dim Mn = 1 and Tr Q = 0, so that from Eq. (1.3) 

1 — e*^ 1 
X(A(.)) = 1- — = 1--. 



Using Eq. (2.8) and letting zq = z — iO+, we get 



H{z) = — Im 
1. 



1 



2^0(^0 - 1/2) 
— Iml 

TT V2O-I/2 Zq 
1 



+ 6{z) 
+ 5{z) 



as expected in Eq. (2.10). 



(2.10) 
(2.11) 

(2.12) 
(2.13) 
(2.14) 



2. Renyi Entropies 



In terms of the spectral density function, the Renyi entropies in Eq. (1.2) are given by 

r-l 



1 

1 - a 







dz ^i(z) ln[2;" + (1 - 2)"]. 



(2.15) 



Dropping the (5-funcion contributions at z = 0, 1 in Eq. (2.8) where the integrand vanishes and 
performing integration by parts by transferring the derivative, we obtain 



Sn 



— ^ / dz- Im[lnx(A(z-iO+))]. 

7r(Q-l)yo 2" + (l-2)"-l ^ ^ ^ ^" 



We can now neglect the small imaginary part 0"*". Let 



Then, after some simplification, 



Sa — 



tanh(Qu) — tanhw , „ , 

au Imllnx(— TT — 2m)\. 

a — 1 



(2.16) 



(2.17) 



(2.18) 



7 



Im z 





ill 
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-ill 



Rez 



FIG. 1: The contour C used in Eq. (2.23). The poles (crosses) are shown for the case o = 5. 



In particular, for the von Neumann entropy 

1 u 

5 = lim 5a = -- / du 2— Im[lnx(-7r - 2m))]. (2.19) 

TT J _^ cosh U 

We can expand the cumulant-generating function as follows: 

lnx(-^ - l.u) = J2 t ^2-2°^ 

k=l 

Now, suppose (contrary to reality) that it is justified to integrate the series term by term, i.e., switch 
the order of summation and integration, as 

Sc, = Y,l3k{a)Ck (2.21) 

k=l 

where 

l3Ja) = — {— Im / du [tanh(Q7ru) - tanh(7rii)] lu . (2.22) 

{l-a)k\ J_^ \ 2 J 

From Im z = (z — z*)/{2i) we see that (3k{o) = for k odd. For k even we can let z = u — i/2 and 
rewrite Eq. (2.22) as the contour integral 

= TT^^^TT^ / [tanh(a7r(z + i/2)) - coth{nz)]z'' (2.23) 
(1 - a)kl 2i Jc 

where C is the rectangular contour shown in Fig. 1. Note that the integrla vanishes on the vertical 
edges where jRe z\ — >• oo. We will compute the integral for a = n,n integer, and analytically continue 
the result to arbitrary a. 

Treating the case of n even and n odd separately, we find the following: 

1. n odd. In this case tanh(n7r(2 + i/2)) = coth(?i7r2) and the integrand has simple poles at 
nz = ±i,±2i, . . . ,±[(n — l)/2]i. The Residue Thereom then gives 

2. n even. In this case tanh{n7r{z + i/2)) — ia,nh{n7rz) and the integrand has simple poles at 
2nz — zbz, dz32, . . . , dz(n — 1). The Residue Theorem gives 

2 1 2mV^ / 1^ 



p=i 
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The two cases can be combined and analytically continued in n to non-integer values: 



2 1 J , Q + 1 



where (^{s,a) is the Hurwitz zeta function, cf. Eq. (4.10). 

It is instructive to check that Eq. (2.26) reduces to the result of Ref. 3 for a = 1. To take the limit 
as a —> 1 we can use I'Hopital's rule and the identities 



daC{s,a) = -s(;{s + l,a), (2.27) 

C(i-s) 

where r(a;) is the gamma function. Then 



C(l -s) = 2(27r)-^r(.s)C(s) cos ^ (2.28) 



hm /?,(«) = lim 2(^a„cf - k, (2.29) 



C(l-fc) (2.30) 

(2.31) 



(fc-1)! 

f2C(fc) if fc even, 
|o if fc odd, 



which agrees with the coefRcients given in Ref. [3], recalling that for even k we have 2({k) = 
(27r)*|i3fc|/fc! where Bk are the Bernoulli numbers. 

Of course, the expansion in Eq. (2.20) is not valid in the context of integrating Eq. (2.18) term 
by term, which we have done, because the range of the integral extends along the entire real line. 
With the exception of gaussian generating functions of the form In x('^) = ~ o"^ A^/2 with only Ci 
and C2 nonzero (i.e., not an infinite series) all cumulant generating functions occuring in fermionic 
systems possess singularities in the complex plane so that the series in Eq. (2.20) has a finite radius of 
convergence. The result is a subtle one: the coefficients obtained from the previous computation are 
the limiting values of the coefficients of the convergent series derived in the main text (with cutoff), 
but using the series in this limiting form will result in divergent values for the Renyi entanglement 
entropies, including the von Neumann entropy. 



III. DIAGONALIZATION OF THE SPIN-1/2 XX CHAIN 

In this section we derive the correlation matrix for the Hamiltonian 



(3.1) 



for both periodic boundary conditions (PBCs) and open boundary conditions (OBCs). The compu- 
tation is very standard, but is provided to make the discussion in the main text more transparent. 
The spin- 1/2 XX chain is solved by using the Jordan- Wigner transformation [4] 



n(l - 2a]a,) 



J<1 



s- = 



n(l-2a]a,: 



to map Eq. (3.1) to the free-fermion Hamiltonian 



XX 



H ■ 



' 1 



Sf = 



I =a,a^ 



1 

1 ) 



(3.2) 



(3.3) 
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For PBCs the cyclic term induces a term of the form 

SlSf + SlSf = -l{a{a, + al&L) exp (itt ^ a]a,) , (3.4) 



and since for the ground state we are interested in the case of zero total magnetization Sj = 0, 
i.e., half-filling of the fermions = L/2, we have 

' for OBCs, 
1 for PBCs, L = 2 mod 4, (3.5) 
-1 for PBCs, L = mod 4. 

We will always assume that L is even. We hope that the dual role of i as an index and as the 
imaginary unit will not cause confusion. 



A. PBCs 

For PBCs the eigenvectors cf) and eigenvalues J A satisfy 

^<^,_i - A0, + ^cf>j+, = 0, (3.6) 

subject to the boundary condition cAl+i = f^i- The eigenvalues and corresponding eigenvectors are 
Afc = cosfc, (t)jk = e'^^V-v/Z where 

, , , 27r n if p = 1. . , 

k(n) = — x{ , / ■ 3.7 

^ ^ L |n + l/2 ifp = -l ^ ' 



for n = —L/2 + 1, . . . , L/2. The correlation matrix (1.8) is therefore 

L/2 

Ml''''{P) = j^ J2 e^'^'-^^fiPJAk). (3.8) 

n=-L/2+l 

Although the sum in Eq. (3.8) cannot be carried out analytically for general L and /3, there are two 
limits of interest where it is possible. First, at zero temperature we have 

/(/3JA,)^''^-.(-A,) = |0 (3.9) 

so that Eq. (3.8) becomes, for both p = ±1, 

fPBC/« 7 . _N _ X 1 sinc[7r(i - j)/2] 



2 smc[7r(i — j)/LJ 



where 



In particular, for L — > oo 



1 if a; = 0, 
otherwise. 



Mf^'^{pJ, L^oo)= S,., - ^sinc 



(3.11) 



(3.12) 
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Second, at high temperature /3J <C 1 

f{liJAk « 1) 
so that the correlation matrix is given by 



l_PJ_ 
2 T 



Afc + 0((/3J)3), 



L/2 

^ e**(^--') cos k 

n=-L/2+l 



+p5\i-jlL-l) 



(3.13) 

(3.14) 
(3.15) 



B. OBCs 

For OBCs the eigenvectors and eigenvalues JA also satisfy Eq. (3.6) but the boundary conditions 
imply = 0L+1 = so that 



At = cos k. 



L+l 



sin jfc, 



where k{n) = ttti/ (L + 1) for n = 1, . . . ,L. The correlation matrix is therefore 



M°«°(/3) = ^ sin(jfc) sin(jfc)/(;3JA,) 

n=l 

1 ^ 

= E cos[(i - j>]/(/3^A,) - {{t - j) + j)}. 



At zero temperature we can use Eq. (3.9) to obtain 

MOI'C(^j ^ oo) = - ^g{i - j, L)sinc 



where the finite-size correction factor is 



8\L+lJ sinc{7rj7[2(L + l)]} 
and goes to 1 — cos(7rj'/2) for L — > co. For L — > oo we therefore have 



+ ^gi-i + j, i)sinc 



1/ TTjL y sinc^{7rjL/[4(L + 1)]} 



M°''''ij3J,L^oo) = S,j 



-M-j) 



+ -smc 



[i + i) 



(3.16) 

(3.17) 
(3.18) 



(3.19) 



(3.20) 



(3.21) 



IV. FLUCTUATION SUMS 

In this section we show how to compute spin fluctuations in one dimension, 



{s-){s-)]. 



(4.1) 
(4.2) 
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In all of the cases considered here the on-site term is 1/4, {{S?)^) — {Sf)'^ = 1/4, so that it will be 
convenient to organize the sums as 

m = i + E (4.3) 

where 

- - . '^' = ^^ (4.4) 
(SfS"!) - (-Sf)(S'p otherwise. ^ ' 

Before presenting the calculations we list some mathematical formulas and notation needed in the 
remainder of the section. 



A. Mathematical Formulas 

We will make extensive use of the digamma function ■^(a;) = d[lnr(x)]/iix, the logarithmic deriva- 
tive of the gamma function r(a;). The digamma function has the series representation 



(4.5) 



and integral representation 



e-* 



V(x)=/ (4.6) 



Some special values of the digamma function and its derivatives are 



V'(l) = -7, M-)=-^-2\n2, ^'(1) = ^- = 

^ ^ k=i 



1 -K^ 



(4.7) 



The digamma function is particularly useful in converting sums to a form that can be expanded using 
the asymptotic expansion 

^(.^oo)=ln.---^^ (4.8) 
= ln^-^-^ + 0(a=-^), (4.9) 

where B„ are the Bernoulli numbers. 

For sums involving non-integer exponents a useful generalization is the Hurwitz zeta function 

oo ^ 

which at a = 1 reduces to the Riemann zeta function 



1 

c(s)=c(s,i)=5:-. (4.11) 

The Htuwitz zeta function has the integral representation 

1 



1 /-oo ,s-y -at 

C(.,a) = ^/ ^i^^. (4.12) 
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The derivative with respect to the second argument is given by 

daC{s,a) = -sC{s + l,a). (4.13) 

We wiU make use of the asymptotic expansion of the Hurwitz zeta function in its second argument 
[5]: 

« - - (73^ + i + 12^ + ^'-''^ 
Double sums that depend only on the difference of the indices can be converted to a single sum as 

e t 

fi\i-j\) = fm + 2Y,{£-k)f{k), (4.15) 

i,j=l k=l 

while sums that depend only on the sum of the indices can be written as 

e e 

J2 fi' + J') = + 1) + - k)f{k + l + e)]. (4.16) 

i,j=l k=l 



fc=i k=i 
In the following computations we will use the notation 

I 1 if fc is odd, 
U II fc IS even. 



B. Spin- 1/2 XX Chain 

For the spin-1/2 XX chain with zero magnetization (Sf) = 0, Wick's Theorem leads to 

- = - |A/,,f , (4.19) 

where the correlation matrix Mij is computed in Sec. III. We therefore carry out Eq. (4.3) with 

Py=-|M,,f, i^j. (4.20) 

1. PBCs 



Here we show how to compute the sum in Eq. (4.3) for [cf. Eq. (3.12)] 



X(7,_,. (4.21) 



Using Eq. (4.15) the required sum is 



13 



where, using Eq. (4.17), the sum on the right-hand side can be rewritten as 



~k (fc + £)2 



OO ^ ^ CXD 



k=l 



k k + ai + (./a 



(4.23) 
(4.24) 



where for the second term in square brackets we have used the trick of converting ak+i to by 
shifting k in the denominator. By writing the sums as (odd terms) = (aU terms) — (even terms) we 
can use Eq. (4.5) to obtain 



A;2 ~ 8 



-V^l 1 + — + 

2 U 2 2q' 



= -^-^-^(l + a,+£) + i0fl + ^ + ^ 
8 2 ^^ ^ -^ 2^12 2 



Expanding for £ 3> 1 using Eq. (4.9) then gives 



£ - |(ln^ + l + 7 + ln2)+0(r^), 



so that 



^'-^ix W = ln^ + /i, /i = l+7 + ln2. 



(4.25) 



(4.26) 



(4.27) 



(4.28) 



plus 0(i'~2) corrections. This result was actually first obtained in Ref. [6] in an attempt to define 
an effective temperature for the subsystem, though without providing the size of the next-order 
correction or the corresponding result for OBCs presented in the next section. In the main text, 
moreover, an alternative method using the theory of ToepHtz matrices is presented for the derivation 
of Eq. (4.28) 



2. OBCs 



Here we show how to compute the sum in Eq. (4.3) for [cf. Eq. (3.21)] 



where 



Prj + Pij + Pr] 



1 1 

TT^ (i — j)2 

1 1 

TT^ (i + _/)2 



'i+V 



-+ _ 2 (-1)' 

P^3 - 7r2(l2-j2) 

Comparing to Eq. (4.21), we can write 



(4.29) 

(4.30) 
(4.31) 
(4.32) 

(4.33) 
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with 



fc=l 



Using Eq. (4.16) we can write as 
and using Eq. (4.17) we get 

k = \ 

-Alt 



(k + i + iy 



■2(fc- 1) 



O'k+l 



(4.34) 



(4.35) 



(4.36) 



k=l 

+ 2 



0-fc 



0-fc 



fc+l/a-l fc + (l + 2£)/Q-l 



A: + cr£ + (1 + £)/q - 1 



(4.37) 



By adding and subtracting a^/k for each term and noting that (odd terms) = (all terms) - 
(even terms), we can build the series representation of the digamma function, Eq. (4.5), to obtain 



Q 



l + 2£\ 1/1 l + 2i 



2^ 2a J 



27/; ai + 



1 + ae 1 + e 

H 

2 2q 



(4.38) 



V(l + 2i) + +i) + 2V'(1 + fJf + £) - -0 1 + 



1 + 21 

(1 + 2^)V''(1 + 2e) + -^^'{l + + 2(1 + ^)V''(1 + ae + e) 
2 ^ V 2 2 



Expanding for ^ ^ 1 using Eq. (4.9) then gives 



7T^ 1+7 



■In^- 



3 1 

4£ 



+ 0{t 



Turning now to tlie cross-term / |_ in Eq. (4.33), we use 



1 / 1 



j2 _ j2 



2i \ j -i j + i 
to write (note the cancellation of the terms with i — j = even) 



E— E(-^-^)+E 



-j=i 



1 



1 



(4.39) 



(4.40) 



(4.41) 



(4.42) 
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The expression in square brackets is 

3 1 

so that we can write 



with 



+ iP{i + e-i)-^{i + £ + i), 



A + B 



2 ^ I- 

' =1 



B = + e-i)-i;{i + e + i) 



(4.43) 
(4.44) 

(4.45) 
(4.46) 



We evaluate each term in turn. For A in Eq. (4.45) we again rewrite the sum as an infinite series, 
Eq. (4.17), then use (even terms) — (odd terms) = 2(even terms) — (all terms) for the alternating 
sum: 



(-ir (-1) 



i+e-\ 



j2 {i + £)2 



I + a 



1 



a=e 

1 



i + a i + a/2 



(4.47) 
(4.48) 
(4.49) 



Adding and subtracing 1/i in the sum to build tlie series representation of the digamma function in 
Eq. (4.5), we obtain 



A- 



8 2^ ' 



1^(1 + 0) -ip\l + 



(4.50) 
(4.51) 



Upon expanding the derivatives of the digamma functions for £ ^ 1 using Eq. (4.9) only the constant 
term contributes to order £~^, so that 



A 



8 



+ 0{£-^). 



(4.52) 



Meanwhile, it has not been possible to give a closed-form expression for B in Eq. (4.46), even in 
terms of digamma functions. We can, however, directly find the asymptotic behavior for large £. To 
this end we first use the integral representation of 't/j{x), Eq. (4.6), to write 



B 



dt" 



1 - e 



i=l 



(-e*)' _ (-e-*)' 
i i 



After writing 



y- = y--y — = -in(i-,T) - ^v — 



i+t 
^l+i/I' 



(4.53) 



(4.54) 
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we expand the factor (1 + i/i) ^ in the second term to get 

2=1 k=0 1=1 

Let 



Then we can write B as 



with 



°° (■— 



k=0 



-{i+e)t i + e* 



In: 



1 + e- 



and 



Ck 



-{i+i)t 

dt ^ ^lfk{-e') - fk{-e-% 

i — e ' 



(4.55) 



(4.56) 



(4.57) 



(4.58) 



(4.59) 



Before evaluating these terms, we note that some of the previous manipulations are only true formally 
because of the range of a; = e* in the integral of Eq. (4.53). However, the purpose is to exclude terms 
proportional to or smaller, which, as we will see, we have accomphshed. 

Now, noting that the argument of the logarithm is simply e*, we have from Eqs. (4.12) and (4.14) 
that 



Next, 



dt 



-{i+i)t 



1 - e 



-t 



,(l+£)t -{l+i)t 



1 + e* 1 + e 



-t 



dt 



-t g-2(l+£)t 



1 - e-2t 1 



dt 

7/'(l + £)-^J 



-t „-(i+£)t 



t 1 - e-* 
2 



dt 



-t/2 



t 1 - e-* 



where we have once again used Eq. (4.6). Expanding for \ with Eq. (4.9), we obtain 
{-if+^CQ = ^[V'(l+^)+7 + 21n2] = ^(ln£ + 7 + 21n2) + 0(ri). 



The next term is 



dt 



g-(i+f)t 
1 - e-* 



g(l+£)t 



-(1+^)4 



(l + e')2 (1 + , 



-t^2 



(4.60) 

(4.61) 
(4.62) 
(4.63) 
(4.64) 

(4.65) 

(4.66) 
(4.67) 
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and gives, after a similar calculation, 

(-l/+ici = Jln£ + 0(r2), (4.68) 
while Cfc for fc > 1 can be shown to only yield terms of 0{£~^) or smaller. Thus 

^5^(0 ^ i^jr (20 + 5^ - ^lln(2q + , + ...21 

plus 0{£~'^) corrections, recalling the result in Eq. (4.28) for PBCs. As one might expect from 
conformal field theory arguments [7] , up to subleading corrections the OBC result is exactly half the 
PBC result with the substitution £ 2£. 



C. Spin-1/2 XXZ Chain (Luttinger Liquids) 



Here we show how to compute the sum in Eq. (4.3) for the Luttinger Liquid result 



^1 , f^„2„,^x-2 4 (-!)'"'• 



27r2 



I* - J 



(4.70) 



We assume that 1/2 < if < 1, which corresponds to anisotropy < A < 1 in the spin-1/2 XXZ 
Hamiltonian where this correlation function is valid. As we will see, the only subleading correction 
relevant to the fluctuations is the term containing Ai. We will, however, keep all of the terms for 
completeness to show that this is indeed the case. Using the formulas derived here it may also be 
possible in the future to compute the constant term and demonstrate explicitly the vanishing of an 
^-dependent term if the non-universal coefficients Am are known [8] . 
Using Eq. (4.15) the sum over Eq. (4.70) becomes 



£ K °° 
-Fxxz(^) = T - -^Foi£) + V 22-'^-i^^F,„(£) 

m=l 



with 



Fo(^) = i:(^-fc)p 



Fm{£) = J2^£-k) 



(4.71) 



(4.72) 



(4.73) 



The usual method of rewriting the finite sums as infinite sums using Eq. (4.17) and converting the 
resulting sums to the digamma function, Eq. (4.5), gives 



CXD OC 



k=i fc=i 

6 Q 



k {k + £f 

f fi ' 

^ k + £/a 



a 



1+ - +7 



£-^(^il + £)-£iP'{l + £)--f. 



(4.74) 
(4.75) 

(4.76) 
(4.77) 
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Expanding for £ 2> 1 using Eq. (4.9), we obtain 



Fo{() = - lii^ - (1 + 7) + 0{r^). 



For -Fi(f) we can similarly write 



^ ( 1 \?7J 



k=l k=l 

mk 



(-1) 



/t = l 



(-ir 



2m'^K - 1 



l-2m^A' \ " 

^(k + 



\ i k + i) 



{-ir 



2m? K 



(4.78) 



(4.79) 



(4.80) 



We must now treat the case of even and odd m separately for each term. In the first term of Eq. (4.80) 
for even m the sum is simply the Riemann zeta function, Eq. (4.11), while for odd m we use the 
identity (even terms) — (odd terms) = 2(even terms) — (all terms) to get 



E 

k=l 



C,(2m?K) for m even, 

^^2i-2m^K _ i^(^(2m^K) for m odd. 



(4.81) 



In the same way, the second term of Eq. (4.80) can also be evaluated as 



2-^ U2m^K-l 



k=l 



C,{2m?K — 1) for m even, 

^^22-2m^K _ i)(^{2m^K - 1) for m odd. 



(4.82) 



In the third term of Eq. (4.80), the sum for m even is the Hurwitz zeta function in Eq. (4.10), 



(-ir 



2m^K- 



- =C 2771^7^-1,1 + 



(?77 even) 



(4.83) 



so that 



1 

277127^-1' 



(k + £/a) 



= -(;{2m^K-l,l+i) 



mk 

^ (fc + £/q)2'""'^-i 

e 



2m?K - 1 



dM 2m^K - 1,1 + - 



= -((277127^ -1,1+1)+ £(;{2m^K, 1 + i) 



(4.84) 
(4.85) 



where in the last line we have used Eq. (4.13) for the derivative of the Hurwitz zeta function. By 
expanding for £ ^ 1 using Eq. (4.14) it can be seen that for even m (smallest 771 = 2 and K > 1/2) 
the sum does not contribute to the overall result up to 0(£~^) corrections. 
Meanwhile, for m odd we simply separate the even and odd terms to write 



(-1)"^ 



^ (/c + f/a)2'»^-^-i 



~,2m''K-l 



C 2m^K-l,l + 



2a 



C(2m^K -1,1 + 



2a 



(m odd). 



(4.86) 
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Then 



1 



(-1) ' 
-(fc + £/Q) 



C( 2m'K - 1, 1 + - ) + C( 2m'K " 1, 3 + 3 



1 



£ 



1 e 



2 2 



(4.87) 



where we have again made use of the derivative formula (4.13). In this case, expanding for £ ^ 1 
using Eq. (4.14) gives corrections smaller than 0{£~^) for all m except m = 1, in which case 



I d 



S (fe + 



1 1 



+ 0{t 



m = 1. 



Combining the previous results, we find 

' (21-2'"' -f^ _ i)c^{2m''K)t - (22-2-'-ff - l)C(2m2A- - 1) 

■^-L m-l 
(2i-2m^if - l)C(2m2i4')£ - (22-2'^'^ _ I)((2m2i4' - 1) m odd, m > 1, 
X(2m2i4')£ - ((2m2_R' - 1) m even. 

+ 0(^2), 



(4.88) 



so that 



^'-Fxxz(^) =/^lnf + 



(-1)^ 

23-2ii: pK 



(4.89) 
(4.90) 



plus 0(£ 2) corrections, where 



/ll 1 + 7 



171=2,4,6,... 



^ (2-22'" ■^-i)C(2m2i4'- 1)^^ 

m=l,3,5,... 

2m'^K-l, 



'C{2m^K - l)A„ 



(4.91) 



and we have dropped a linear term oll^ on the assumption that 

,2m=x-i)^(2m2/i:)yl„ 



fflLL 



m=l,3, 
2m^if-l 



C{2m^K)A„ 



0. 



(4.92) 
(4.93) 



If all of the coefficients Am are known this assumption can, in principle, be demonstrated explicitly. 
We can check, for example, that the first few terms are correct at the free-fermion point A = where 



1 

2^' 



An 



0, 



(4.94) 



which follows from comparing Eq. (4.70) with Eq. (4.21). Then in Eq. (4.92) clearly oll = 0, while 
in Eq. (4.91) 



III = 1+1 -I lim(2 - 22-i)C(2x - 1) = 1 + 7 + In 2, 



(4.95) 



consistent with Eq. (4.28). In the next section a more elaborate example with slightly different 
exponents is given. 
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D. Haldane-Shastry Chain 



Here we show how to compute the sum in Eq. (4.3) for [9-11] 

1 , Si(7rr) , , , ^ 

Pv = 4(-ir^, r-\i-j\- (4.96) 

Before doing an exact calculation we may observe that pij has the asymptotic expansion 

1 1 (-l)"- 1 ^ (-l)™(2m)! 



lyizinM! (4 97) 



Comparison to Eq. (4.70) shows that the leading behavior of the fluctuations can be obtained by 
letting 



in Eq. (4.90), which gives 



K=l, A, = ^ (4.98) 



^'-^Hs(^) = ^ln£ + /HS-^4^ (4.99) 



plus 0{£~'^) corrections. The value K = 1/2 is expected from the fact that the Haldane-Shastry 
model is in the same universality class as the isotropic Heisenberg chain, for which /•i' = 1/2 at the 
renormalization group flxed point. The difference is that the correlation functions of the Haldane- 
Shastry model do not have logarithmic corrections due to marginal (in the renormalization group 
sense) operators. 

For the Haldane-Shastry model we can also show explicitly that the "volume" term aus^ corre- 
sponding to Eq. (4.92) vanishes, and moreover compute the numerical value of the constant term 
corresponding to Eq. (4.91). Taking the difference in exponents between Eq. (4.70) and Eq. (4.96) 
into account, we find 

aHs_l-ln2 1 ^ (-l)"'+i(2m- 2)! 



r2 4 2 

m=l 



("00) 

and 

iE^=i5l^a.-i.. (4.101, 



l^ = ^ + 1±2 + (-l)'"(2m)!^ 
7r2 8 27r2 27r2 

ni=l 



Both sums appear to be divergent since the Riemann zeta function ({x) — > 1 as a; — > co so that 
terms grow larger in magnitude: this is not surprising since we summed an asymptotic expansion. 
Nevertheless, we can show that the following divergent summation methods produce the correct 
expressions for the coefficient of the linear term (zero) and constant term, as compared to numerical 
results from DMRG presented in the main text. 

First, recalling that (-l)™+i(2m)!C(2m)/7r2'" = 22'"-i52m where B„ are the Bernoulli numbers, 
we can write Eq. (4.100) as 

OHs ^ l-ln2 1 B^^ 

7r2 4 2 (2m)(2m - 1) (l/2)2™-i ^ 



l-ln2 
4 2 



^ CO 

E 



771=1 



, i'2777 

ax 



1/2 



(2m)a;2™ 



(4.103) 
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Reversing the sum and integral and comparing to the series representation of the digamma function, 
Eq. (4.8), we see that 



OHs ^ l-hi2 _ 1 



dx 



1/2 



In a; ^(x) 

2.T 



(4.104) 



Using the original definition of the digamma function, tlj{x) = d[\nT{x)]/dx, and lnr(a;) — > {x 
1/2) In a; — X + (1/2) ln(27r) for x — > oo, we arrive at the expected result 



ohs 1 — In 2 1 
^ ~ 4 2 



X ■ 



In a; — a; — lnr(a;) 



1/2 



l-ln2 l/l-ln2 
4 2 

0. 



Meanwhile, the series in Eq. (4.101) can be carried out using Borel summation. Let 

oo 



k=l 



with 



The desired value is then 



The Borel transform of y{z) is 



i^k\C^[k + 1) for k even, 
~ for k odd. 



/hs ^ 1 l + i + y{l/iT) 
7r2 8 27r2 



^2'W = E!f*' = E^(2fc + l)(it) = 



(4.105) 

(4.106) 
(4.107) 



(4.108) 



(4.109) 



(4.110) 



(4.111) 



fc=i 



There is a useful identity, essentially the Taylor expansion of the digamma function at z = 1, that 
relates the digamma function with the Riemann zeta function: 



^C(fc + i)(-^)' = -7-V'(i + ^). 



k=l 



The even terms can be singled out by using the alternating factor (— 1)'^: 



Bv(t) = -7 - -[^{l + it) + V(l - it)]. 



The Borel sum of y{z) is thus 



y{z) 



dt e-*[By{tz) 



7 / dt e"* 



1 



dt [ip{l + itz) + V'(l - itz)] 



dt [tp(l + itz) + - itz)]. 



(4.112) 

(4.113) 

(4.114) 
(4.115) 
(4.116) 
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We can now write Eq. (4.110) as 



/hs ^ 1 I- A 
7r2 8 27r2 



where 



A-- 



dt e 



(4.117) 



(4.118) 



It is worth rewriting the constant A to keep all expressions real. Using the integral representation of 
tp{x), Eq. (4.G), and switching the order of integration, we obtain 



A: 



dx 



dx e 



dt e 



3 ^ e ^ xt 

^ cos — 

X 1 — e ^ TT 



X 1 — e ^ 1 + (x/7r)2 



Numerical integration of Eq. (4.120) yields 



/] 



HS 



0.197217. 



(4.119) 
(4.120) 

(4.121) 



V. THE RELATION BETWEEN FACTORIAL AND ORDINARY CUMULANTS 

In this section we derive the general relation between factorial and ordinary cumulants, which is 
not new but to our knowledge is not well-known. The ordinary cumulants, or simply "cumulants," of 
a probability distribution P{X) are defined as the coefficients C,i in the expansion of the cumulant 
generating function 



where xl-*^) is the characteristic function 

x(A) = (e^^^) = ^e'^^P(X). 

X 

Thus 

Cn = {-idxrinxW\x=o- 



(5.1) 



(5.2) 



(5.3) 



Meanwhile, the factorial cumulants of the same probability distribution are defined as the coeflacients 
Fn in the expansion of the factorial cumulant generating function 



n=l 

where Xf{^) is the factorial moment generating function 

X/(^) = ((^ + l)^) = ^(^ + l)^P(X). 

X 

Therefore 



(5.4) 



(5.5) 



d:inxf{z)U=o. 



(5.6) 
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We see that the two generating functions xW and Xfi^) a-re related by 2; + 1 = e*'^, i.e., 

X/W=xHln(^+l)), (5.7) 
x(A)=x/(e^'-l). (5.8) 

The two relations (5.7) and (5.8) can be used to express the two kinds of cumulants in terms of 
each other. First, to express the factorial cumulants in terms of ordinary cumulants we use Eq. (5.7) 
and the series representations (5.1) and (5.4) to obtain 



k=l 



k\ 



The RHS of Eq. (5.9) can be rewritten using the identity [12] 

^ s{n,k) „ ^ [ln(^ + 1)]^ 



(5.9) 



(5.10) 



where s(n, k) are the signed Stirhng numbers of the first kind, so that 



s{n, k) 



n=l k=l Ti=k 

Switching the order of the sums on the RHS gives 



00 ^ r n 



n=l n=l '- 

and by equating coefficients we obtain 

n 

F„ = y~]s(n,fc)Cfc. 



k=l 



(5.11) 



(5.12) 



(5.13) 



Similarly, to express the ordinary cumulants in terms of factorial cumulants we first write, using 
Eq. (5.8); 



= 1^ — ki — 

11=1 k=l 



The identity [12] 



^ S{n,k) ^^^^,,_{e^^-lf 

71= k 



k\ 



(5.14) 



(5.15) 



where S{n, k) are the Stirling numbers of the second kind allows us to rewrite the RHS of Eq. (5.14) 
to get 



(5.16) 



and switching the order of sums gives 



-k=l 



(5.17) 
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Equating coefficients, we then obtain 



Cn = Y,S{n,k)Fk. (5.18) 
fc=i 

Eq. (5.18) can also be obtained from Eq. (5.13) by noting that s{n,k) and S{n,k) are inverses of 
each other when considered as square matrices. 

For convenience, in the main text we use the notation Si{n, k) = |s(n, k)\ = (— l)"~'^.s(n, k) for the 
unsigned Stirhng numbers of the first kind. 
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